Markov state models of biomolecular conformational dynamics

Self-organizing maps of unbiased ligand-target binding pathways and kinetics

The interpretation of ligand-target interactions at atomistic resolution is central to most efforts in computational drug discovery and optimization. However, the highly dynamic nature of protein targets, as well as possible induced fit effects, makes difficult to sample many interactions effectively with docking studies or even with large-scale molecular dynamics (MD) simulations. We propose a novel application of Self-Organizing Maps (SOMs) to address the sampling and dynamic mapping tasks, particularly in cases involving ligand flexibility and induced fit. The SOM approach offers a data-driven strategy to create a map of the interaction process and pathways based on unbiased MD. Furthermore, we show how the preliminary SOM mapping is complementary to kinetic analysis, with the employment of both network-based approaches and Markov state models. We demonstrate the method by comprehensively mapping a large dataset of 640 μs of unbiased trajectories sampling the recognition process between the phosphorylated YEEI peptide and its high-specificity target lck-SH2. The integration of SOM into unbiased simulation protocols significantly advances our understanding of the ligand binding mechanism. This approach serves as a potent tool for mapping intricate ligand-target interactions with unprecedented detail, thereby enhancing the characterization of kinetic properties crucial to drug design.

CDR L3 Loop Rearrangement Switches Multispecific SPE-7 IgE Antibody From Hapten to Protein Binding

The monoclonal IgE antibody SPE-7 was originally raised against a 2,4-dinitrophenyl (DNP) target. Through its ability to adopt multiple conformations, the antibody is capable of binding to a diverse range of small haptens and large proteins. The present study examines a dataset of experimentally determined crystal structures of the SPE-7 antibody to gain insight into the mechanisms that contribute to its multispecificity. With the emergence of more and more therapeutic antibodies against a huge repertoire of different targets, our research could be of great interest for future drug development. We are able to discriminate between the different paratope-binding states in the conformational ensembles obtained by enhanced sampling molecular dynamics simulations, and to calculate their transition timescales and state probabilities. Furthermore, we describe the key residues responsible for discriminating between the different binding capacities and identify a tryptophan in a central position of the CDR L3 loop as the residue of greatest interest. The overall dynamics of the paratope appear to be mainly influenced by the CDR L3 and CDR L1 loops.

Protein representations: Encoding biological information for machine learning in biocatalysis

Enzymes offer a more environmentally friendly and low-impact solution to conventional chemistry, but they often require additional engineering for their application in industrial settings, an endeavour that is challenging and laborious. To address this issue, the power of machine learning can be harnessed to produce predictive models that enable the in silico study and engineering of improved enzymatic properties. Such machine learning models, however, require the conversion of the complex biological information to a numerical input, also called protein representations. These inputs demand special attention to ensure the training of accurate and precise models, and, in this review, we therefore examine the critical step of encoding protein information to numeric representations for use in machine learning. We selected the most important approaches for encoding the three distinct biological protein representations - primary sequence, 3D structure, and dynamics - to explore their requirements for employment and inductive biases. Combined representations of proteins and substrates are also introduced as emergent tools in biocatalysis. We propose the division of fixed representations, a collection of rule-based encoding strategies, and learned representations extracted from the latent spaces of large neural networks. To select the most suitable protein representation, we propose two main factors to consider. The first one is the model setup, which is influenced by the size of the training dataset and the choice of architecture. The second factor is the model objectives such as consideration about the assayed property, the difference between wild-type models and mutant predictors, and requirements for explainability. This review is aimed at serving as a source of information and guidance for properly representing enzymes in future machine learning models for biocatalysis.

Assessing AF2's ability to predict structural ensembles of proteins

Recent breakthroughs in protein structure prediction have enhanced the precision and speed at which protein configurations can be determined. Additionally, molecular dynamics (MD) simulations serve as a crucial tool for capturing the conformational space of proteins, providing valuable insights into their structural fluctuations. However, the scope of MD simulations is often limited by the accessible timescales and the computational resources available, posing challenges to comprehensively exploring protein behaviors. Recently emerging approaches have focused on expanding the capability of AlphaFold2 (AF2) to predict conformational substates of protein. Here, we benchmark the performance of various workflows that have adapted AF2 for ensemble prediction and compare the obtained structures with ensembles obtained from MD simulations and NMR. We provide an overview of the levels of performance and accessible timescales that can currently be achieved with machine learning (ML) based ensemble generation. Significant minima of the free energy surfaces remain undetected.

AI in cellular engineering and reprogramming

During the last decade, artificial intelligence (AI) has increasingly been applied in biophysics and related fields, including cellular engineering and reprogramming, offering novel approaches to understand, manipulate, and control cellular function. The potential of AI lies in its ability to analyze complex datasets and generate predictive models. AI algorithms can process large amounts of data from single-cell genomics and multiomic technologies, allowing researchers to gain mechanistic insights into the control of cell identity and function. By integrating and interpreting these complex datasets, AI can help identify key molecular events and regulatory pathways involved in cellular reprogramming. This knowledge can inform the design of precision engineering strategies, such as the development of new transcription factor and signaling molecule cocktails, to manipulate cell identity and drive authentic cell fate across lineage boundaries. Furthermore, when used in combination with computational methods, AI can accelerate and improve the analysis and understanding of the intricate relationships between genes, proteins, and cellular processes. In this review article, we explore the current state of AI applications in biophysics with a specific focus on cellular engineering and reprogramming. Then, we showcase a couple of recent applications where we combined machine learning with experimental and computational techniques. Finally, we briefly discuss the challenges and prospects of AI in cellular engineering and reprogramming, emphasizing the potential of these technologies to revolutionize our ability to engineer cells for a variety of applications, from disease modeling and drug discovery to regenerative medicine and biomanufacturing.

Engineering salt-tolerant Cas12f1 variants for gene-editing applications

CRISPR has revolutionized the field of genome editing in life sciences by serving as a versatile and state-of-the-art tool. Cas12f1 is a small nuclease of the bacterial immunity CRISPR system with an ideal size for cellular delivery, in contrast to CRISPR-associated (Cas) proteins like Cas9 or Cas12. However, Cas12f1 works best at low salt concentrations. In this study, we find that the plasticity of certain Cas12f1 regions (K196-Y202 and I452-L515) is negatively affected by increased salt concentrations. On this line, key protein domains (REC1, WED, Nuc, lid) that are involved in the DNA-target recognition and the activation of the catalytic RuvC domain are in turn also affected. We suggest that salt concentration should be taken in to consideration for activity assessments of Cas engineered variants, especially if the mutations are on the protospacer adjacent motif interacting domain. The results can be exploited for the engineering of Cas variants and the assessment of their activity at varying salt concentrations. We propose that the K198Q mutation can restore at great degree the compromised plasticity and could potentially lead to salt-tolerant Cas12f1 variants. The methodology can be also employed for the study of biomolecules in terms of their salinity tolerance.Communicated by Ramaswamy H. Sarma.

Insight into structural dynamics involved in activation mechanism of full length KRAS wild type and P-loop mutants

KRAS protein is known to be frequently mutated in various cancers. The most common mutations being at position 12, 13 and 61. The positions 12 and 13 form part of the phosphate binding region (P-loop) of KRAS. Owing to mutation, the protein remains in continuous active state and affects the normal cellular process. Understanding the structural changes owing to mutations in GDP-bound (inactive state) and GTP-bound (active state) may help in the design of better therapeutics. To understand the structural flexibility due to the mutations specifically located at P-loop regions (G12D, G12V and G13D), extensive molecular dynamics simulations (24 μs) have been carried for both inactive (GDP-bound) and active (GTP-bound) structures for the wild type and these mutants. The study revealed that the local structural changes at the site of mutations allosterically guide changes in distant regions of the protein through hydrogen bond and hydrophobic signalling network. The dynamic cross correlation analysis and the comparison of the correlated motions among different systems manifested that changes in SW-I, SW-II, α3 and the loop preceding α3 affects the interactions of GDP/GTP with different regions of the protein thereby affecting its hydrolysis. Further, the Markov state modelling analysis confirmed that the mutations, especially G13D imparts rigidity to structure compared to wild type and thus limiting its conformational state in either intermediate state or active state. The study suggests that along with SW-I and SW-II regions, the loop region preceding the α3 helix and α3 helix are also involved in affecting the hydrolysis of nucleotides and may be considered while designing therapeutics against KRAS.

The effect of the loop on the thermodynamic and kinetic of single base pair in pseudoknot

RNA pseudoknots are RNA molecules with specialized three-dimensional structures that play important roles in various biological processes. To understand the functions and mechanisms of pseudoknots, it is essential to elucidate their structures and folding pathways. The most fundamental step in RNA folding is the opening and closing of a base pair. The effect of flexible loops on the base pair in pseudoknots remains unclear. In this work, we use molecular dynamics simulations and Markov state model to study the configurations, thermodynamic and kinetic of single base pair in pseudoknots. We find that the presence of the loop leads to a trap state. In addition, the rate-limiting step for the formation of base pair is the disruption of the trap state, rather than the open state to the closed state, which is quite different from the previous studies on non-pseudoknot RNA. For the thermodynamic parameters in pseudoknots, we find that the entropy difference upon opening the base pair between this simulation and the nearest-neighbor model results from the different entropy of different lengths of loop in solution. The thermodynamic parameters of the stack in pseudoknot are close to the nearest-neighbor parameters. The bases on the loop have different distribution patterns in different states, and the slow transition states of the loop are determined by the orientation of the bases.

The fast committor machine: Interpretable prediction with kernels

In the study of stochastic systems, the committor function describes the probability that a system starting from an initial configuration x will reach a set B before a set A. This paper introduces an efficient and interpretable algorithm for approximating the committor, called the "fast committor machine" (FCM). The FCM uses simulated trajectory data to build a kernel-based model of the committor. The kernel function is constructed to emphasize low-dimensional subspaces that optimally describe the A to B transitions. The coefficients in the kernel model are determined using randomized linear algebra, leading to a runtime that scales linearly with the number of data points. In numerical experiments involving a triple-well potential and alanine dipeptide, the FCM yields higher accuracy and trains more quickly than a neural network with the same number of parameters. The FCM is also more interpretable than the neural net.

Molecular Dynamics Simulation Combined with Neural Relationship Inference and Markov Model to Reveal the Relationship between Conformational Regulation and Bioluminescence Properties of Gaussia Luciferase

Gaussia luciferase (Gluc) is currently known as the smallest naturally secreted luciferase. Due to its small molecular size, high sensitivity, short half-life, and high secretion efficiency, it has become an ideal reporter gene and is widely used in monitoring promoter activity, studying protein-protein interactions, protein localization, high-throughput drug screening, and real-time monitoring of tumor occurrence and development. Although studies have shown that different Gluc mutations exhibit different bioluminescent properties, their mechanisms have not been further investigated. The purpose of this study is to reveal the relationship between the conformational changes of Gluc mutants and their bioluminescent properties through molecular dynamics simulation combined with neural relationship inference (NRI) and Markov models. Our results indicate that, after binding to the luciferin coelenterazine (CTZ), the α-helices of the 109-119 residues of the Gluc Mutant2 (GlucM2, the flash-type mutant) are partially unraveled, while the α-helices of the same part of the Gluc Mutant1 (GlucM1, the glow-type mutant) are clearly formed. The results of Markov flux analysis indicate that the conformational differences between glow-type and flash-type mutants when combined with luciferin substrate CTZ mainly involve the helicity change of α7. The most representative conformation and active pocket distance analysis indicate that compared to the flash-type mutant GlucM2, the glow-type mutant GlucM1 has a higher degree of active site closure and tighter binding. In summary, we provide a theoretical basis for exploring the relationship between the conformational changes of Gluc mutants and their bioluminescent properties, which can serve as a reference for the modification and evolution of luciferases.

Modulation of the conformational landscape of the PDZ3 domain by perturbation on a distal non-canonical α3 helix: decoding the microscopic mechanism of allostery in the PDZ3 domain

While allosteric signal transduction is crucial for protein signaling and regulation, the dynamic process of allosteric communication remains poorly understood. The third PDZ domain (PDZ stands for the common structural domain shared by the postsynaptic density protein (PSD95), Drosophila disc large tumor suppressor (DlgA), and zonula occludens-1 protein (ZO-1)) serves as a classic example of a single-domain allosteric protein, demonstrating a long-range coupling between the C-terminal α helix (known as the α3 helix) and ligand binding. A molecular level understanding of how the α3 helix modulates the ligand binding affinity of the PDZ3 domain is still lacking. In this study, extensive molecular dynamics simulations corroborated with principal component analysis (PCA), ligand binding free energy calculations, energetic frustration analysis and Markov state model analysis are employed to uncover such molecular details. We demonstrate the definite presence of a binding competent closed-like state in the conformational landscape of wild-type PDZ3. The population modulations of this closed state and other binding incompetent states in the landscape due to α3-truncation/mutation of PDZ3 are explored. A correlation between the closed state population and calculated binding free energy is established, which supports the conformation selection mechanism. Covariance analysis identified the presence of correlated motion between two distant loops (β1-β2 and β2-β3) in the wild-type PDZ3 system, which weakened due to truncation/mutation in the distant α3 helix. It has also been observed that whenever the α3 helix was perturbed, the β2-β3 loop got further away from the binding groove and it is found to be correlated with the binding free energy values. Energetic frustration analysis of the PDZ3 domain also showed that the β2-β3 loop is highly frustrated. Finally, MSM analysis revealed a relevant timescale (closed to open state transition), which is similar to the observed experimental signal transduction timescale for the system. These observations led to the conclusion that the distantly located α3 helix plays a pivotal role in regulating the conformational landscape of the PDZ3 domain, determining the ligand binding affinity and resulting in allosteric behavior of the domain.

Pathways and mechanism of MRTX1133 binding to KRAS G12D elucidated by molecular dynamics simulations and Markov state models

KRAS G12D is the most common oncogenic mutation identified in several types of cancer. Therefore, design of inhibitors targeting KRAS G12D represents a promising strategy for anticancer therapy. MRTX1133 is a highly potent inhibitor (approximate experiment Kd ≈ 0.0002 nM) of KRAS G12D and is currently in Phase 1/2 study, however, pathways of the compound binding to KRAS G12D has remained unknown, and the mechanism underlying the complicated dynamic process are challenging to capture experimentally, which hinder the structure-based anti-cancer drug design. Here, using MRTX1133 as a probe, unbiased molecular dynamics (MD) was used to simulate the process of MRTX1133 spontaneously binding to KRAS G12D. In six of 42 independent MD simulation (a total of 99 μs), MRTX1133 was observed to successfully associate with KRAS G12D. The kinetically metastable states refer to the potential pathways of MRTX1133 binding to KRAS G12D were revealed by Markov state models (MSM) analysis. Additionally, 8 key residues that are essential for MRTX1133 recognition and tight binding at the preferred low energy states were identified by MM/GBSA analysis. In sum, this study provides a new perspective on understanding the pathways and mechanism of MRTX1133 binding to KRAS G12D.

Interactions between Inhibitors and 5-Lipoxygenase: Insights from Gaussian Accelerated Molecular Dynamics and Markov State Models

Inflammation is a protective stress response triggered by external stimuli, with 5-lipoxygenase (5LOX) playing a pivotal role as a potent mediator of the leukotriene (Lts) inflammatory pathway. Nordihydroguaiaretic acid (NDGA) functions as a natural orthosteric inhibitor of 5LOX, while 3-acetyl-11-keto-β-boswellic acid (AKBA) acts as a natural allosteric inhibitor targeting 5LOX. However, the precise mechanisms of inhibition have remained unclear. In this study, Gaussian accelerated molecular dynamics (GaMD) simulation was employed to elucidate the inhibitory mechanisms of NDGA and AKBA on 5LOX. It was found that the orthosteric inhibitor NDGA was tightly bound in the protein's active pocket, occupying the active site and inhibiting the catalytic activity of the 5LOX enzyme through competitive inhibition. The binding of the allosteric inhibitor AKBA induced significant changes at the distal active site, leading to a conformational shift of residues 168-173 from a loop to an α-helix and significant negative correlated motions between residues 285-290 and 375-400, reducing the distance between these segments. In the simulation, the volume of the active cavity in the stable conformation of the protein was reduced, hindering the substrate's entry into the active cavity and, thereby, inhibiting protein activity through allosteric effects. Ultimately, Markov state models (MSM) were used to identify and classify the metastable states of proteins, revealing the transition times between different conformational states. In summary, this study provides theoretical insights into the inhibition mechanisms of 5LOX by AKBA and NDGA, offering new perspectives for the development of novel inhibitors specifically targeting 5LOX, with potential implications for anti-inflammatory drug development.

Mutual Linearity of Nonequilibrium Network Currents

For continuous-time Markov chains and open unimolecular chemical reaction networks, we prove that any two stationary currents are linearly related upon perturbations of a single edge's transition rates, arbitrarily far from equilibrium. We extend the result to nonstationary currents in the frequency domain, provide and discuss an explicit expression for the current-current susceptibility in terms of the network topology, and discuss possible generalizations. In practical scenarios, the mutual linearity relation has predictive power and can be used as a tool for inference or model proof testing.

Unveiling Allosteric Regulation and Binding Mechanism of BRD9 through Molecular Dynamics Simulations and Markov Modeling

Bromodomain-containing protein 9 (BRD9) is a key player in chromatin remodeling and gene expression regulation, and it is closely associated with the development of various diseases, including cancers. Recent studies have indicated that inhibition of BRD9 may have potential value in the treatment of certain cancers. Molecular dynamics (MD) simulations, Markov modeling and principal component analysis were performed to investigate the binding mechanisms of allosteric inhibitor POJ and orthosteric inhibitor 82I to BRD9 and its allosteric regulation. Our results indicate that binding of these two types of inhibitors induces significant structural changes in the protein, particularly in the formation and dissolution of α-helical regions. Markov flux analysis reveals notable changes occurring in the α-helicity near the ZA loop during the inhibitor binding process. Calculations of binding free energies reveal that the cooperation of orthosteric and allosteric inhibitors affects binding ability of inhibitors to BRD9 and modifies the active sites of orthosteric and allosteric positions. This research is expected to provide new insights into the inhibitory mechanism of 82I and POJ on BRD9 and offers a theoretical foundation for development of cancer treatment strategies targeting BRD9.

Capturing the Effects of Single Atom Substitutions on the Inhibition Efficiency of Glycogen Synthase Kinase-3β Inhibitors via Markov State Modeling and Experiments

Small modifications in the chemical structure of ligands are known to dramatically change their ability to inhibit the activity of a protein. Unraveling the mechanisms that govern these dramatic changes requires scrutinizing the dynamics of protein-ligand binding and unbinding at the atomic level. As an exemplary case, we have studied Glycogen Synthase Kinase-3β (GSK-3β), a multifunctional kinase that has been implicated in a host of pathological processes. As such, there is a keen interest in identifying ligands that inhibit GSK-3β activity. One family of compounds that are highly selective and potent inhibitors of GSK-3β is exemplified by a molecule termed COB-187. COB-187 consists of a five-member heterocyclic ring with a thione at C2, a pyridine substituted methyl at N3, and a hydroxyl and phenyl at C4. We have studied the inhibition of GSK-3β by COB-187-related ligands that differ in a single heavy atom from each other (either in the location of nitrogen in their pyridine ring, or with the pyridine ring replaced by a phenyl ring), or in the length of the alkyl group joining the pyridine and the N3. The inhibition experiments show a large range of half-maximal inhibitory concentration (IC50) values from 10 nM to 10 μM, implying that these ligands exhibit vastly different propensities to inhibit GSK-3β. To explain these differences, we perform Markov State Modeling (MSM) using fully atomistic simulations. Our MSM results are in excellent agreement with the experiments in that they accurately capture differences in the binding propensities of the ligands. The simulations show that the binding propensities are related to the ligands' ability to attain a compact conformation where their two aromatic rings are spatially close. We rationalize this result by sampling numerous binding and unbinding events via funnel metadynamics simulations, which show that indeed while approaching the bound state, the ligands prefer to be in their compact conformation. We find that the presence of nitrogen in the aromatic ring increases the probability of attaining the compact conformation. Protein-ligand binding is understood to be dictated by the energetics of interactions and entropic factors, like the release of bound water from the binding pockets. This work shows that changes in the conformational distribution of ligands due to atom-level modifications in the structure play an important role in protein-ligand binding.

Identification of a Cryptic Pocket in Methionine Aminopeptidase-II Using Adaptive Bandit Molecular Dynamics Simulations and Markov State Models

Methionine aminopeptidase-II (MetAP-II) is a metalloprotease, primarily responsible for the cotranslational removal of the N-terminal initiator methionine from the nascent polypeptide chain during protein synthesis. MetAP-II has been implicated in angiogenesis and endothelial cell proliferation and is therefore considered a validated target for cancer therapeutics. However, there is no effective drug available against MetAP-II. In this study, we employ Adaptive Bandit molecular dynamics simulations to investigate the structural dynamics of the apo and ligand-bound MetAP-II. Our results focus on the dynamic behavior of the disordered loop that is not resolved in most of the crystal structures. Further analysis of the conformational flexibility of the disordered loop reveals a hidden cryptic pocket that is predicted to be potentially druggable. The network analysis indicates that the disordered loop region has a direct signaling route to the active site. These findings highlight a new way to target MetAP-II by designing inhibitors for the allosteric site within this disordered loop region.

Using Gaussian accelerated molecular dynamics combined with Markov state models to explore the mechanism of action of new oral inhibitors on Complex I

In this study, our focus was on investigating H-1,2,3-triazole derivative HP661 as a novel and highly efficient oral OXPHOS inhibitor, with its molecular-level inhibitory mechanism not yet fully understood. We selected the ND1, NDUFS2, and NDUFS7 subunits of Mitochondrial Complex I as the receptor proteins and established three systems for comparative analysis: protein-IACS-010759, protein-lead compound 10, and protein-HP661. Through extensive analysis involving 500 ns Gaussian molecular dynamics simulations, we gained insights into these systems. Additionally, we constructed a Markov State Models to examine changes in secondary structures during the motion processes. The research findings suggest that the inhibitor HP661 enhances the extensibility and hydrophilicity of the receptor protein. Furthermore, HP661 induces the unwinding of the α-helical structure in the region of residues 726-730. Notably, key roles were identified for Met37, Phe53, and Pro212 in the binding of various inhibitors. In conclusion, we delved into the potential molecular mechanisms of triazole derivative HP661 in inhibiting Complex I. These research outcomes provide crucial information for a deeper understanding of the mechanisms underlying OXPHOS inhibition, offering valuable theoretical support for drug development and disease treatment design.

Interplay between phosphorylation and oligomerization tunes the conformational ensemble of SWEET transporters

SWEET sugar transporters are desirable biotechnological targets for improving plant growth. One engineering strategy includes modulating how SWEET transporters are regulated. Phosphorylation and oligomerization have been shown to positively regulate SWEET function, leading to increased sugar transport activity. However, constitutive phosphorylation may not be beneficial to plant health under basal conditions. Structural and mechanistic understanding of the interplay between phosphorylation and oligomerization in functional regulation of SWEETs remains limited. Using extensive molecular dynamics simulations coupled with Markov state models, we demonstrate the thermodynamic and kinetic effects of SWEET phosphorylation and oligomerization using OsSWEET2b as a model. We report that the beneficial effects of these SWEET regulatory mechanisms bias outward-facing states and improved extracellular gating, which complement published experimental findings. Our results offer molecular insights to SWEET regulation and may guide engineering strategies throughout the SWEET transport family.

Secondary Structure Stabilization of Macrocyclic Antimicrobial Peptides via Cross-Link Swapping

Lactam cross-links have been employed to stabilize the helical secondary structure and enhance the activity and physiological stability of antimicrobial peptides; however, stabilization of β-sheets via lactamization has not been observed. In the present study, lactams between the side chains of C- and N-terminal residues have been used to stabilize the β-sheet conformation in a short ten-residue analogue of chicken angiogenin-4. Designed using a combination of molecular dynamics simulations and Markov state models, the lactam cross-linked peptides are shown to adopt stabilized β-sheet conformations consistent with simulated structures. Replacement of the peptide side-chain Cys-Cys disulfide by a lactam cross-link enhanced the broad-spectrum antibacterial activity compared to the parent peptide and exhibited greater propensity to induce proinflammatory activity in macrophages. The combination of molecular simulations and conformational and biological analyses of the synthetic peptides provides a useful paradigm for the rational design of therapeutically active peptides with constrained β-sheet structures.

Testing the Simplified Molecular Dynamics Approach to Improve the Reproduction of ECD Spectra and Monitor Aggregation

A simplified molecular-dynamics-based electronic circular dichroism (ECD) approach was tested on three condensed derivatives with limited conformational flexibility and an isochroman-2H-chromene hybrid, the ECD spectra of which could not be precisely reproduced by the conventional ECD calculation protocol. Application of explicit solvent molecules at the molecular mechanics (MD) level in the dynamics simulations and subsequent TDDFT-ECD calculation for the unoptimized MD structures was able to improve the agreements between experimental and computed spectra. Since enhancements were achieved even for molecules with limited conformational flexibility, deformations caused by the solvent molecules and multitudes of conformers produced with unoptimized geometries seem to be key factors for better agreement. The MD approach could confirm that aggregation of the phenanthrene natural product luzulin A had a significant contribution to a specific wavelength range of the experimental ECD. The MD approach has proved that dimer formation occurred in solution and this was responsible for the anomalous ECD spectrum. The scope and limitations of the method have also been discussed.

Information thermodynamics of transition paths between multiple mesostates

A central concern across the natural sciences is a quantitative understanding of the mechanism governing rare transitions between two metastable states. Recent research has uncovered a fundamental equality between the time-reversal asymmetry of the ensemble of such transition paths and the informativeness of system dynamics about the reactivity of a given trajectory, immediately leading to quantitative criteria for judging the importance of distinct system coordinates for the transition. Here we generalize this framework to multiple mesostates. We find that the main system-wide and coordinate-specific results generalize intuitively, while the combinatorial diversity of pairwise transitions raises new questions and points to new concepts. This work increases the previous framework's generality and applicability and forges connections to enhanced-sampling and coarse-grained dynamical approaches such as milestoning and Markov-state models.

Conformational Dynamics in Corynebacterium glutamicum Diaminopimelate Epimerase: Insights from Ligand Parameterization, Atomistic Simulation, and Markov State Modeling

The microbial enzyme diaminopimelate epimerase (DapF), a vital enzyme in the lysine biosynthetic pathway, catalyzes the conversion of L, L-diaminopimelate (L, L-DAP) to D, L-diaminopimelate (D, L-DAP) using a catalytic cysteine dyad with one cysteine in thiol state and another in thiolate. Under oxidizing conditions, the catalytic cysteines of apo DapF form a disulfide bond that alters the structure and function of DapF. Given its potential as a target for antimicrobial resistance treatments, understanding DapF's functional dynamics is imperative. In the present work, we employ microsecond-scale all-atom molecular dynamics simulations of product-bound DapF and apo-DapF under oxidized and reduced conditions. We employ a polarized charge model for the ligand and the active site residues, which was necessary to preserve the electrostatic environment in the active site and retain the ligand in the active site. The product-bound DapF and apo-DapF in oxidized and reduced conditions exhibit a closed, semi-open, and open conformation, respectively, as identified using the internal coordinates of the dimeric enzyme and the principal component analysis. The conformational switch is guided by the dynamic catalytic (DC) loop, loop II, and loop III movements in the active site. The time scale of the close-to-open conformational transition is estimated to be 0.8 μs through Markov state modeling (MSM) and transition path theory (TPT). The present study explains the role of various active site residues and loops in ligand binding and protein dynamics in the DapF enzyme under different redox conditions. Such information will be helpful in future inhibitor design studies targeting the DapF enzyme.

Quantitative Prediction of Dissociation Rates of PYK2 Ligands Using Umbrella Sampling and Milestoning

We used umbrella sampling and the milestoning simulation method to study the dissociation of multiple ligands from protein kinase PYK2. The activation barriers obtained from the potential of mean force of the umbrella sampling simulations correlated well with the experimental dissociation rates. Using the zero-temperature string method, we obtained optimized paths along the free-energy surfaces for milestoning simulations of three ligands with a similar structure. The milestoning simulations gave an absolute dissociation rate within 2 orders of magnitude of the experimental value for two ligands but at least 3 orders of magnitude too high for the third. Despite the similarity in their structures, the ligands took different pathways to exit from the binding site of PYK2, making contact with different sets of residues. In addition, the protein experienced different conformational changes for dissociation of the three ligands.

A kinetic model reveals the critical gating motifs for donor-substrate loading into Actinobacillus pleuropneumoniae N-glycosyltransferase

Soluble N-glycosyltransferase from Actinobacillus pleuropneumoniae (ApNGT) catalyzes the glycosylation of asparagine residues, and represents one of the most encouraging biocatalysts for N-glycoprotein production. Since the sugar tolerance of ApNGT is restricted to limited monosaccharides (e.g., Glc, GlcN, Gal, Xyl, and Man), tremendous efforts are devoted to expanding the substrate scope of ApNGT via enzyme engineering. However, rational design of novel NGT variants suffers from an elusive understanding of the substrate-binding process from a dynamic point of view. Here, by employing extensive all-atom molecular dynamics (MD) simulations integrated with a kinetic model, we reveal, at the atomic level, the complete donor-substrate binding process from the bulk solvent to the ApNGT active-site, and the key intermediate states of UDP-Glc during its loading dynamics. We are able to determine the critical transition event that limits the overall binding rate, which guides us to pinpoint the key ApNGT residues dictating the donor-substrate entry. The functional roles of several identified gating residues were evaluated through site-directed mutagenesis and enzymatic assays. Two single-point mutations, N471A and S496A, could profoundly enhance the catalytic activity of ApNGT. Our work provides deep mechanistic insights into the structural dynamics of the donor-substrate loading process for ApNGT, which sets a rational basis for design of novel NGT variants with desired substrate specificity.

Understanding the Fast-Triggering Unfolding Dynamics of FK-11 upon Photoexcitation of Azobenzene

Photoswitchable molecules can control the activity and functions of biomolecules by triggering conformational changes. However, it is still challenging to fully understand such fast-triggering conformational evolution from nonequilibrium to equilibrium distribution at the molecular level. Herein, we successfully simulated the unfolding of the FK-11 peptide upon the photoinduced trans-to-cis isomerization of azobenzene based on the Markov state model. We found that the ensemble of FK-11 contains five conformational states, constituting two unfolding pathways. More intriguingly, we observed the microsecond-scale conformational propagation of the FK-11 peptide from the fully folded state to the equilibrium populations of the five states. The computed CD spectra match well with the experimental data, validating our simulation method. Overall, our study not only offers a protocol to study the photoisomerization-induced conformational changes of enzymes but also could orientate the rational design of a photoswitchable molecule to manipulate biological functions.

Unraveling motion in proteins by combining NMR relaxometry and molecular dynamics simulations: A case study on ubiquitin

Nuclear magnetic resonance (NMR) relaxation experiments shine light onto the dynamics of molecular systems in the picosecond to millisecond timescales. As these methods cannot provide an atomically resolved view of the motion of atoms, functional groups, or domains giving rise to such signals, relaxation techniques have been combined with molecular dynamics (MD) simulations to obtain mechanistic descriptions and gain insights into the functional role of side chain or domain motion. In this work, we present a comparison of five computational methods that permit the joint analysis of MD simulations and NMR relaxation experiments. We discuss their relative strengths and areas of applicability and demonstrate how they may be utilized to interpret the dynamics in MD simulations with the small protein ubiquitin as a test system. We focus on the aliphatic side chains given the rigidity of the backbone of this protein. We find encouraging agreement between experiment, Markov state models built in the χ1/χ2 rotamer space of isoleucine residues, explicit rotamer jump models, and a decomposition of the motion using ROMANCE. These methods allow us to ascribe the dynamics to specific rotamer jumps. Simulations with eight different combinations of force field and water model highlight how the different metrics may be employed to pinpoint force field deficiencies. Furthermore, the presented comparison offers a perspective on the utility of NMR relaxation to serve as validation data for the prediction of kinetics by state-of-the-art biomolecular force fields.

The Role of Force Fields and Water Models in Protein Folding and Unfolding Dynamics

Protein folding is a fascinating, not fully understood phenomenon in biology. Molecular dynamics (MD) simulations are an invaluable tool to study conformational changes in atomistic detail, including folding and unfolding processes of proteins. However, the accuracy of the conformational ensembles derived from MD simulations inevitably relies on the quality of the underlying force field in combination with the respective water model. Here, we investigate protein folding, unfolding, and misfolding of fast-folding proteins by examining different force fields with their recommended water models, i.e., ff14SB with the TIP3P model and ff19SB with the OPC model. To this end, we generated long conventional MD simulations highlighting the perks and pitfalls of these setups. Using Markov state models, we defined kinetically independent conformational substates and emphasized their distinct characteristics, as well as their corresponding state probabilities. Surprisingly, we found substantial differences in thermodynamics and kinetics of protein folding, depending on the combination of the protein force field and water model, originating primarily from the different water models. These results emphasize the importance of carefully choosing the force field and the respective water model as they determine the accuracy of the observed dynamics of folding events. Thus, the findings support the hypothesis that the water model is at least equally important as the force field and hence needs to be considered in future studies investigating protein dynamics and folding in all areas of biophysics.

Main role of fractal-like nature of conformational space in subdiffusion in proteins

Protein dynamics involves a myriad of mechanical movements happening at different time and space scales, which make it highly complex. One of the less understood features of protein dynamics is subdiffusivity, defined as sublinear dependence between displacement and time. Here, we use all-atoms molecular dynamics (MD) simulations to directly interrogate an already well-established theory and demonstrate that subdiffusivity arises from the fractal nature of the network of metastable conformations over which the dynamics, thought of as a diffusion process, takes place.

Toward physics-based precision medicine: Exploiting protein dynamics to design new therapeutics and interpret variants

The goal of precision medicine is to utilize our knowledge of the molecular causes of disease to better diagnose and treat patients. However, there is a substantial mismatch between the small number of food and drug administration (FDA)-approved drugs and annotated coding variants compared to the needs of precision medicine. This review introduces the concept of physics-based precision medicine, a scalable framework that promises to improve our understanding of sequence-function relationships and accelerate drug discovery. We show that accounting for the ensemble of structures a protein adopts in solution with computer simulations overcomes many of the limitations imposed by assuming a single protein structure. We highlight studies of protein dynamics and recent methods for the analysis of structural ensembles. These studies demonstrate that differences in conformational distributions predict functional differences within protein families and between variants. Thanks to new computational tools that are providing unprecedented access to protein structural ensembles, this insight may enable accurate predictions of variant pathogenicity for entire libraries of variants. We further show that explicitly accounting for protein ensembles, with methods like alchemical free energy calculations or docking to Markov state models, can uncover novel lead compounds. To conclude, we demonstrate that cryptic pockets, or cavities absent in experimental structures, provide an avenue to target proteins that are currently considered undruggable. Taken together, our review provides a roadmap for the field of protein science to accelerate precision medicine.

Functional regulation of aquaporin dynamics by lipid bilayer composition

With the diversity of lipid-protein interactions, any observed membrane protein dynamics or functions directly depend on the lipid bilayer selection. However, the implications of lipid bilayer choice are seldom considered unless characteristic lipid-protein interactions have been previously reported. Using molecular dynamics simulation, we characterize the effects of membrane embedding on plant aquaporin SoPIP2;1, which has no reported high-affinity lipid interactions. The regulatory impacts of a realistic lipid bilayer, and nine different homogeneous bilayers, on varying SoPIP2;1 dynamics are examined. We demonstrate that SoPIP2;1's structure, thermodynamics, kinetics, and water transport are altered as a function of each membrane construct's ensemble properties. Notably, the realistic bilayer provides stabilization of non-functional SoPIP2;1 metastable states. Hydrophobic mismatch and lipid order parameter calculations further explain how lipid ensemble properties manipulate SoPIP2;1 behavior. Our results illustrate the importance of careful bilayer selection when studying membrane proteins. To this end, we advise cautionary measures when performing membrane protein molecular dynamics simulations.

Diffusive dynamics of a model protein chain in solution

A Markov state model is a powerful tool that can be used to track the evolution of populations of configurations in an atomistic representation of a protein. For a coarse-grained linear chain model with discontinuous interactions, the transition rates among states that appear in the Markov model when the monomer dynamics is diffusive can be determined by computing the relative entropy of states and their mean first passage times, quantities that are unchanged by the specification of the energies of the relevant states. In this paper, we verify the folding dynamics described by a diffusive linear chain model of the crambin protein in three distinct solvent systems, each differing in complexity: a hard-sphere solvent, a solvent undergoing multi-particle collision dynamics, and an implicit solvent model. The predicted transition rates among configurations agree quantitatively with those observed in explicit molecular dynamics simulations for all three solvent models. These results suggest that the local monomer-monomer interactions provide sufficient friction for the monomer dynamics to be diffusive on timescales relevant to changes in conformation. Factors such as structural ordering and dynamic hydrodynamic effects appear to have minimal influence on transition rates within the studied solvent densities.

Quantifying Disorder in a Protein by Mapping its Locally Correlated Structure and Kinetics

Proteins under physiological conditions are inherently mobile and sample a vast array of structures. Consequently, the need arises, on the one hand, at a local level to determine the independent moving parts and their associated conformations and kinetics, and on the other hand, at a global level, to quantify the disorder in the full protein molecule. We present an approach that provides these quantities in the form of local kinetic network models, which are constructed by analyzing the molecular dynamics (MD) trajectories of the protein molecule. Entropies of independent parts of the molecule are quantified. The method outlined here, using the Trp-cage miniprotein prototype, offers a new tool to understand the dynamic structural changes that ultimately govern the functioning of a protein. The method is particularly suited to problems where there are subtle changes in the structure or dynamics at local levels, for example, due to ligand binding.

Integrating cellular electron microscopy with multimodal data to explore biology across space and time

Biology spans a continuum of length and time scales. Individual experimental methods only glimpse discrete pieces of this spectrum but can be combined to construct a more holistic view. In this Review, we detail the latest advancements in volume electron microscopy (vEM) and cryo-electron tomography (cryo-ET), which together can visualize biological complexity across scales from the organization of cells in large tissues to the molecular details inside native cellular environments. In addition, we discuss emerging methodologies for integrating three-dimensional electron microscopy (3DEM) imaging with multimodal data, including fluorescence microscopy, mass spectrometry, single-particle analysis, and AI-based structure prediction. This multifaceted approach fills gaps in the biological continuum, providing functional context, spatial organization, molecular identity, and native interactions. We conclude with a perspective on incorporating diverse data into computational simulations that further bridge and extend length scales while integrating the dimension of time.

Broadly inhibitory antibodies against severe malaria virulence proteins

Plasmodium falciparum pathology is driven by the accumulation of parasite-infected erythrocytes in microvessels. This process is mediated by the parasite's polymorphic erythrocyte membrane protein 1 (PfEMP1) adhesion proteins. A subset of PfEMP1 variants that bind human endothelial protein C receptor (EPCR) through their CIDRα1 domains is responsible for severe malaria pathogenesis. A longstanding question is whether individual antibodies can recognize the large repertoire of circulating PfEMP1 variants. Here, we describe two broadly reactive and binding-inhibitory human monoclonal antibodies against CIDRα1. The antibodies isolated from two different individuals exhibited a similar and consistent EPCR-binding inhibition of 34 CIDRα1 domains, representing five of the six subclasses of CIDRα1. Both antibodies inhibited EPCR binding of both recombinant full-length and native PfEMP1 proteins as well as parasite sequestration in bioengineered 3D brain microvessels under physiologically relevant flow conditions. Structural analyses of the two antibodies in complex with two different CIDRα1 antigen variants reveal similar binding mechanisms that depend on interactions with three highly conserved amino acid residues of the EPCR-binding site in CIDRα1. These broadly reactive antibodies likely represent a common mechanism of acquired immunity to severe malaria and offer novel insights for the design of a vaccine or treatment targeting severe malaria.

Sequence - dynamics - function relationships in protein tyrosine phosphatases

Protein tyrosine phosphatases (PTPs) are crucial regulators of cellular signaling. Their activity is regulated by the motion of a conserved loop, the WPD-loop, from a catalytically inactive open to a catalytically active closed conformation. WPD-loop motion optimally positions a catalytically critical residue into the active site, and is directly linked to the turnover number of these enzymes. Crystal structures of chimeric PTPs constructed by grafting parts of the WPD-loop sequence of PTP1B onto the scaffold of YopH showed WPD-loops in a wide-open conformation never previously observed in either parent enzyme. This wide-open conformation has, however, been observed upon binding of small molecule inhibitors to other PTPs, suggesting the potential of targeting it for drug discovery efforts. Here, we have performed simulations of both enzymes and show that there are negligible energetic differences in the chemical step of catalysis, but significant differences in the dynamical properties of the WPD-loop. Detailed interaction network analysis provides insight into the molecular basis for this population shift to a wide-open conformation. Taken together, our study provides insight into the links between loop dynamics and chemistry in these YopH variants specifically, and how WPD-loop dynamic can be engineered through modification of the internal protein interaction network.

Unraveling the allosteric inhibition mechanism of PARP-1 CAT and the D766/770A mutation effects via Gaussian accelerated molecular dynamics and Markov state model

PARP-1 (Poly (ADP-ribose) polymerase 1) is a nuclear enzyme and plays a key role in many cellular functions, such as DNA repair, modulation of chromatin structure, and recombination. Developing the PARP-1 inhibitors has emerged as an effective therapeutic strategy for a growing list of cancers. The catalytic structural domain (CAT) of PARP-1 upon binding the inhibitor allosterically regulates the conformational changes of helix domain (HD), affecting its identification with the damaged DNA. The typical type I (EB47) and III (veliparib) inhibitors were able to lengthening or shortening the retention time of this enzyme on DNA damage and thus regulating the cytotoxicity. Nonetheless, the basis underlying allosteric inhibition is unclear, which limits the development of novel PARP-1 inhibitors. Here, to investigate the distinct allosteric changes of EB47 and veliparib against PARP-1 CAT, each complex was simulated via classical and Gaussian accelerated molecular dynamics (cMD and GaMD). To study the reverse allosteric basis and mutation effects, the complexes PARP-1 with UKTT15 and PARP-1 D766/770A mutant with EB47 were also simulated. Importantly, the markov state models were built to identify the transition pathways of crucial substates of allosteric communication and the induction basis of PARP-1 reverse allostery. The conformational change differences of PARP-1 CAT regulated by allosteric inhibitors were concerned with to their interaction at the active site. Energy calculations suggested the energy advantage of EB47 in inhibiting the wild-type PARP-1, compared with D766/770A PARP-1. Secondary structure results showed the change of two key loops (αB-αD and αE-αF) in different systems. This work reported the basis of PARP-1 allostery from both thermodynamic and kinetic views, providing the guidance for the discovery and design of more innovative PARP-1 allosteric inhibitors.

Revealing the conformational dynamics of UDP-GlcNAc recognition by O-GlcNAc transferase via Markov state model

The O-linked N-acetylglucosamine (O-GlcNAc) glycosylation is a critical post-translational modification and closely linked to various physiological and pathological conditions. The O-GlcNAc transferase (OGT) functions as the only glycosyltransferase of O-GlcNAc glycosylation by transferring GlcNAc from UDP-GlcNAc to serine or threonine residues on protein substrates. The interaction mode of UDP-GlcNAc against OGT has been preliminarily revealed by the crystal structures, yet an atomic-level comprehension for the conformational dynamics of the recognition process remains elusive. Here, we construct the Markov state model based on extensive all-atom molecular dynamics (MD) simulations with an aggregated simulation time of ∼9 μs, and reveal that the UDP-GlcNAc recognition process by OGT encompasses four key metastable states, occurring within an estimated timescale of ∼10 μs. During UDP-GlcNAc recognition process, we find the pyrophosphate moiety (P2O52-) initially anchors to the active pocket via salt bridge and hydrogen bonds, facilitating subsequent binding of the uridine and GlcNAc moieties. Furthermore, the functional roles of K842 involved in the salt bridge with P2O52- were evaluated through extra mutant MD simulations. Overall, our study provides valuable insights into the UDP-GlcNAc recognition mechanism by OGT, which could further aid in mechanistic studies of O-GlcNAc glycosylation and drug development targeting on OGT.

pH-Dependent Structure and Dynamics of the Catalytic Domains of Human Carbonic Anhydrase II and IX

Extensive computer simulation studies have been carried out to probe the pH-dependent structure and dynamics of the two most efficient isoenzymes II and IX of human carbonic anhydrase (HCA) that control the pH in the human body. The equilibrium structure and hydration of their catalytic domains are found to be largely unaffected by the variation of pH in the range studied, in close agreement with the known experimental results. In contrast, a significant effect of the change in pH is observed for the first time on the local electrostatic potential of the active site walls and the dynamics of active site water molecules. We also report for the first time the free energy and kinetics of coupled fluctuations of orientation and protonation states of the well-known His-mediated proton shuttle (His-64) in both isozymes at pH 7 and 8. The transitions between different tautomers of in or out conformations of His-64 side chain range between 109 and 106 s-1 depending on pH. Possible implications of these results on conformation-dependent pKa of His-64 side chain and its role in driving the catalysis toward hydration of CO2 or dehydration of HCO3- with varying pH are discussed.

Phage display assisted discovery of a pH-dependent anti-α-cobratoxin antibody from a natural variable domain library

Recycling IgG antibodies bind to their target antigen at physiological pH in the blood stream and release them upon endocytosis when pH levels drop, allowing the IgG antibodies to be recycled into circulation via FcRn-mediated cellular pathways, while the antigens undergo lysosomal degradation. This enables recycling antibodies to achieve comparable therapeutic effect at lower doses than their non-recycling counterparts. The development of such antibodies is typically achieved by histidine doping of their variable regions or by performing in vitro antibody selection campaigns utilizing histidine doped libraries. Both are strategies that may introduce sequence liabilities. Here, we present a methodology that employs a naïve antibody phage display library, consisting of natural variable domains, to discover antibodies that bind α-cobratoxin from the venom of Naja kaouthia in a pH-dependent manner. As a result, an antibody was discovered that exhibits a 7-fold higher off-rate at pH 5.5 than pH 7.4 in bio-layer interferometry experiments. Interestingly, no histidine residues were found in its variable domains, and in addition, the antibody showed pH-dependent binding to a histidine-devoid antigen mutant. As such, the results demonstrate that pH-dependent antigen-antibody binding may not always be driven by histidine residues. By employing molecular dynamics simulations, different protonation states of titratable residues were found, which potentially could be responsible for the observed pH-dependent antigen binding properties of the antibody. Finally, given the typically high diversity of naïve antibody libraries, the methodology presented here can likely be applied to discover recycling antibodies against different targets ab initio without the need for histidine doping.

QMrebind: incorporating quantum mechanical force field reparameterization at the ligand binding site for improved drug-target kinetics through milestoning simulations

Understanding the interaction of ligands with biomolecules is an integral component of drug discovery and development. Challenges for computing thermodynamic and kinetic quantities for pharmaceutically relevant receptor-ligand complexes include the size and flexibility of the ligands, large-scale conformational rearrangements of the receptor, accurate force field parameters, simulation efficiency, and sufficient sampling associated with rare events. Our recently developed multiscale milestoning simulation approach, SEEKR2 (Simulation Enabled Estimation of Kinetic Rates v.2), has demonstrated success in predicting unbinding (koff) kinetics by employing molecular dynamics (MD) simulations in regions closer to the binding site. The MD region is further subdivided into smaller Voronoi tessellations to improve the simulation efficiency and parallelization. To date, all MD simulations are run using general molecular mechanics (MM) force fields. The accuracy of calculations can be further improved by incorporating quantum mechanical (QM) methods into generating system-specific force fields through reparameterizing ligand partial charges in the bound state. The force field reparameterization process modifies the potential energy landscape of the bimolecular complex, enabling a more accurate representation of the intermolecular interactions and polarization effects at the bound state. We present QMrebind (Quantum Mechanical force field reparameterization at the receptor-ligand binding site), an ORCA-based software that facilitates reparameterizing the potential energy function within the phase space representing the bound state in a receptor-ligand complex. With SEEKR2 koff estimates and experimentally determined kinetic rates, we compare and interpret the receptor-ligand unbinding kinetics obtained using the newly reparameterized force fields for model host-guest systems and HSP90-inhibitor complexes. This method provides an opportunity to achieve higher accuracy in predicting receptor-ligand koff rate constants.

How Robust Is the Ligand Binding Transition State?

For many drug targets, it has been shown that the kinetics of drug binding (e.g., on rate and off rate) is more predictive of drug efficacy than thermodynamic quantities alone. This motivates the development of predictive computational models that can be used to optimize compounds on the basis of their kinetics. The structural details underpinning these computational models are found not only in the bound state but also in the short-lived ligand binding transition states. Although transition states cannot be directly observed experimentally due to their extremely short lifetimes, recent successes have demonstrated that modeling the ligand binding transition state is possible with the help of enhanced sampling molecular dynamics methods. Previously, we generated unbinding paths for an inhibitor of soluble epoxide hydrolase (sEH) with a residence time of 11 min. Here, we computationally modeled unbinding events with the weighted ensemble method REVO (resampling of ensembles by variation optimization) for five additional inhibitors of sEH with residence times ranging from 14.25 to 31.75 min, with average prediction accuracy within an order of magnitude. The unbinding ensembles are analyzed in detail, focusing on features of the ligand binding transition state ensembles (TSEs). We find that ligands with similar bound poses can show significant differences in their ligand binding TSEs, in terms of their spatial distribution and protein-ligand interactions. However, we also find similarities across the TSEs when examining more general features such as ligand degrees of freedom. Together these findings show significant challenges for rational, kinetics-based drug design.

SARS-CoV-2 RNA-Dependent RNA Polymerase Follows Asynchronous Translocation Pathway for Viral Transcription and Replication

Translocation is one essential step for the SARS-CoV-2 RNA-dependent RNA polymerase (RdRp) to exert viral replication and transcription. Although cryo-EM structures of SARS-CoV-2 RdRp are available, the molecular mechanisms of dynamic translocation remain elusive. Herein, we constructed a Markov state model based on extensive molecular dynamics simulations to elucidate the translocation dynamics of the SARS-CoV-2 RdRp. We identified two intermediates that pinpoint the rate-limiting step of translocation and characterize the asynchronous movement of the template-primer duplex. The 3'-terminal nucleotide in the primer strand lags behind due to the uneven distribution of protein-RNA interactions, while the translocation of the template strand is delayed by the hurdle residue K500. Even so, the two strands share the same "ratchet" to stabilize the polymerase in the post-translocation state, suggesting a Brownian-ratchet model. Overall, our study provides intriguing insights into SARS-CoV-2 replication and transcription, which would open a new avenue for drug discoveries.

Dynamics of organophosphate guest encapsulation in heteroleptic coordination cages

Heteroleptic coordination cages allow the design of different host structures that can bind guest molecules within their cavities. In a previous work, the energetics of organophosphate encapsulation in palladium(II)-based heteroleptic coordination cages that differ in terms of their ability to form hydrogen bonds have been investigated [Platzek et al., Endohedrally Functionalized Heteroleptic Coordination Cages for Phosphate Ester Binding, Angew. Chem., Int. Ed. 2022, 61, e2022093]. The present work focuses on the dynamics of this system. Dynamic information is obtained through the application of a Markov state model (MSM) to unbiased multi-microsecond atomistic molecular dynamics simulations of guest binding and release. The MSM reveals that both the bound state and the binding/unbinding pathways are highly dynamic, with different types of interactions mediating the binding of the diphenylphosphate guest. Thus, the simulations highlight the dynamic nature of the nanoconfinement in the host-guest systems, with possible implications for the use of such coordination cages as catalysts.

Identification of Aggregation Mechanism of Acetylated PHF6* and PHF6 Tau Peptides Based on Molecular Dynamics Simulations and Markov State Modeling

The microtubule-associated protein tau (MAPT) has a critical role in the development and preservation of the nervous system. However, tau's dysfunction and accumulation in the human brain can lead to several neurodegenerative diseases, such as Alzheimer's disease, Down's syndrome, and frontotemporal dementia. The microtubule binding (MTB) domain plays a significant, important role in determining the tau's pathophysiology, as the core of paired helical filaments PHF6* (275VQIINK280) and PHF6 (306VQIVYK311) of R2 and R3 repeat units, respectively, are formed in this region, which promotes tau aggregation. Post-translational modifications, and in particular lysine acetylation at K280 of PHF6* and K311 of PHF6, have been previously established to promote tau misfolding and aggregation. However, the exact aggregation mechanism is not known. In this study, we established an atomic-level nucleation-extension mechanism of the separated aggregation of acetylated PHF6* and PHF6 hexapeptides, respectively, of tau. We show that the acetylation of the lysine residues promotes the formation of β-sheet enriched high-ordered oligomers. The Markov state model analysis of ac-PHF6* and ac-PHF6 aggregation revealed the formation of an antiparallel dimer nucleus which could be extended from both sides in a parallel manner to form mixed-oriented and high-ordered oligomers. Our study describes the detailed mechanism for acetylation-driven tau aggregation, which provides valuable insights into the effect of post-translation modification in altering the pathophysiology of tau hexapeptides.

Improving Estimation of the Koopman Operator with Kolmogorov-Smirnov Indicator Functions

It has become common to perform kinetic analysis using approximate Koopman operators that transform high-dimensional timeseries of observables into ranked dynamical modes. The key to the practical success of the approach is the identification of a set of observables that form a good basis on which to expand the slow relaxation modes. Good observables are, however, difficult to identify a priori and suboptimal choices can lead to significant underestimations of characteristic time scales. Leveraging the representation of slow dynamics in terms of Hidden Markov Models (HMM), we propose a simple and computationally efficient clustering procedure to infer surrogate observables that form a good basis for slow modes. We apply the approach to an analytically solvable model system as well as on three protein systems of different complexities. We consistently demonstrate that the inferred indicator functions can significantly improve the estimation of the leading eigenvalues of Koopman operators and correctly identify key states and transition time scales of stochastic systems, even when good observables are not known a priori.

Identifying Selectivity Filters in Protein Biosensor for Ligand Screening

Specialized sensing mechanisms in bacteria enable the identification of cognate ligands with remarkable selectivity in highly xenobiotic-polluted environments where these ligands are utilized as energy sources. Here, via integrating all-atom computer simulation, biochemical assay, and isothermal titration calorimetry measurements, we determine the molecular basis of MopR, a phenol biosensor's complex selection process of ligand entry. Our results reveal a set of strategically placed selectivity filters along the ligand entry pathway of MopR. These filters act as checkpoints, screening diverse aromatic ligands at the protein surface based on their chemical features and sizes. Ligands meeting specific criteria are allowed to enter the sensing site in an orientation-dependent manner. Sequence and structural analyses demonstrate the conservation of this ligand entry mechanism across the sensor class, with individual amino acids along the selectivity filter path playing a critical role in ligand selection. Together, this investigation highlights the importance of interactions with the ligand entry pathway, in addition to interactions within the binding pocket, in achieving ligand selectivity in biological sensing. The findings enhance our understanding of ligand selectivity in bacterial phenol biosensors and provide insights for rational expansion of the biosensor repertoire, particularly for the biotechnologically relevant class of aromatic pollutants.

Functionalized Fullerene Potentially Inhibits SARS-CoV-2 Infection by Modulating Spike Protein Conformational Changes

The disease of SARS-CoV-2 has caused considerable morbidity and mortality globally. Spike proteins on the surface of SARS-CoV-2 allow it to bind with human cells, leading to infection. Fullerenes and their derivatives are promising SARS-CoV-2 inhibitors and drug-delivery vehicles. In this study, Gaussian accelerated molecular dynamics simulations and the Markov state model were employed to delve into the inhibitory mechanism of Fullerene-linear-polyglycerol-b-amine sulfate (F-LGPS) on spike proteins. During the study, it was discovered that fullerene derivatives can operate at the interface of the receptor-binding domain (RBD) and the N-terminal domain (NTD), keeping structural domains in a downward conformation. It was also observed that F-LGPS demonstrated superior inhibitory effects on the XBB variant in comparison to the wild-type variant. This study yielded invaluable insights for the potential development of efficient therapeutics targeting the spike protein of SARS-CoV-2.

Toward a structural identification of metastable molecular conformations

Interpreting high-dimensional data from molecular dynamics simulations is a persistent challenge. In this paper, we show that for a small peptide, deca-alanine, metastable states can be identified through a neural net based on structural information alone. While processing molecular dynamics data, dimensionality reduction is a necessary step that projects high-dimensional data onto a low-dimensional representation that, ideally, captures the conformational changes in the underlying data. Conventional methods make use of the temporal information contained in trajectories generated through integrating the equations of motion, which forgoes more efficient sampling schemes. We demonstrate that EncoderMap, an autoencoder architecture with an additional distance metric, can find a suitable low-dimensional representation to identify long-lived molecular conformations using exclusively structural information. For deca-alanine, which exhibits several helix-forming pathways, we show that this approach allows us to combine simulations with different biasing forces and yields representations comparable in quality to other established methods. Our results contribute to computational strategies for the rapid automatic exploration of the configuration space of peptides and proteins.

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.