Mapping the Substrate Recognition Pathway in Cytochrome P450

Kinetic Ensemble of Tau Protein through the Markov State Model and Deep Learning Analysis

The ordered assembly of Tau protein into filaments characterizes Alzheimer's and other neurodegenerative diseases, and thus, stabilization of Tau protein is a promising avenue for tauopathies therapy. To dissect the underlying aggregation mechanisms on Tau, we employ a set of molecular simulations and the Markov state model to determine the kinetics of ensemble of K18. K18 is the microtubule-binding domain of Tau protein and plays a vital role in the microtubule assembly, recycling processes, and amyloid fibril formation. Here, we efficiently explore the conformation of K18 with about 150 μs lifetimes in silico. Our results observe that all four repeat regions (R1-R4) are very dynamic, featuring frequent conformational conversion and lacking stable conformations, and the R2 region is more flexible than the R1, R3, and R4 regions. Additionally, it is worth noting that residues 300-310 in R2-R3 and residues 319-336 in R3 tend to form sheet structures, indicating that K18 has a broader functional role than individual repeat monomers. Finally, the simulations combined with Markov state models and deep learning reveal 5 key conformational states along the transition pathway and provide the information on the microsecond time scale interstate transition rates. Overall, this study offers significant insights into the molecular mechanism of Tau pathological aggregation and develops novel strategies for both securing tauopathies and advancing drug discovery.

Coenzyme-binding pathway on glutamate dehydrogenase suggested from multiple-binding sites visualized by cryo-electron microscopy

The structure of hexameric glutamate dehydrogenase (GDH) in the presence of the coenzyme nicotinamide adenine dinucleotide phosphate (NADP) was visualized using cryogenic transmission electron microscopy to investigate the ligand-binding pathways to the active site of the enzyme. Each subunit of GDH comprises one hexamer-forming core domain and one nucleotide-binding domain (NAD domain), which spontaneously opens and closes the active-site cleft situated between the two domains. In the presence of NADP, the potential map of GDH hexamer, assuming D3 symmetry, was determined at a resolution of 2.4 Å, but the NAD domain was blurred due to the conformational variety. After focused classification with respect to the NAD domain, the potential maps interpreted as NADP molecules appeared at five different sites in the active-site cleft. The subunits associated with NADP molecules were close to one of the four metastable conformations in the unliganded state. Three of the five binding sites suggested a pathway of NADP molecules to approach the active-site cleft for initiating the enzymatic reaction. The other two binding modes may rarely appear in the presence of glutamate, as demonstrated by the reaction kinetics. Based on the visualized structures and the results from the enzymatic kinetics, we discussed the binding modes of NADP to GDH in the absence and presence of glutamate.

3site Multisubstrate-Bound State of Cytochrome P450cam

We identified a multisubstrate-bound state, hereby referred as a 3site state, in cytochrome P450cam via integrating molecular dynamics simulation with nuclear magnetic resonance (NMR) pseudocontact shift measurements. The 3site state is a result of simultaneous binding of three camphor molecules in three locations around P450cam: (a) in a well-established "catalytic" site near heme, (b) in a kink-separated "waiting" site along channel-1, and (c) in a previously reported "allosteric" site at E, F, G, and H helical junctions. These three spatially distinct binding modes in the 3site state mutually communicate with each other via homotropic allostery and act cooperatively to render P450cam functional. The 3site state shows a significantly superior fit with NMR pseudo contact shift (PCS) data with a Q-score of 0.045 than previously known bound states and consists of D251 free of salt-bridges with K178 and R186, rendering the enzyme functionally primed. To date, none of the reported cocomplex of P450cam with its redox partner putidaredoxin (pdx) has been able to match solution NMR data and controversial pdx-induced opening of P450cam's channel-1 remains a matter of recurrent discourse. In this regard, inclusion of pdx to the 3site state is able to perfectly fit the NMR PCS measurement with a Q-score of 0.08 and disfavors the pdx-induced opening of channel-1, reconciling previously unexplained remarkably fast hydroxylation kinetics with a koff of 10.2 s-1. Together, our findings hint that previous experimental observations may have inadvertently captured the 3site state as an in vitro solution state, instead of the catalytic state alone, and provided a distinct departure from the conventional understanding of cytochrome P450.

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.

Reversing the Enantioselectivity of Enzymatic Carbene N-H Insertion Through Mechanism-Guided Protein Engineering

We report a computationally driven approach to access enantiodivergent enzymatic carbene N-H insertions catalyzed by P411 enzymes. Computational modeling was employed to rationally guide engineering efforts to control the accessible conformations of a key lactone-carbene (LAC) intermediate in the enzyme active site by installing a new H-bond anchoring point. This H-bonding interaction controls the relative orientation of the reactive carbene intermediate, orienting it for an enantioselective N-nucleophilic attack by the amine substrate. By combining MD simulations and site-saturation mutagenesis and screening targeted to only two key residues, we were able to reverse the stereoselectivity of previously engineered S-selective P411 enzymes. The resulting variant, L5_FL-B3, accepts a broad scope of amine substrates for N-H insertion with excellent yields (up to >99 %), high efficiency (up to 12 300 TTN), and good enantiocontrol (up to 7 : 93 er).

Markov State Models Reconcile Conformational Plasticity of GTPase with Its Substrate Binding Event

Ras GTPase is an enzyme that catalyzes the hydrolysis of guanosine triphosphate (GTP) and plays an important role in controlling crucial cellular signaling pathways. However, this enzyme has always been believed to be undruggable due to its strong binding affinity with its native substrate GTP. To understand the potential origin of high GTPase/GTP recognition, here we reconstruct the complete process of GTP binding to Ras GTPase via building Markov state models (MSMs) using a 0.1 ms long all-atom molecular dynamics (MD) simulation. The kinetic network model, derived from the MSM, identifies multiple pathways of GTP en route to its binding pocket. While the substrate stalls onto a set of non-native metastable GTPase/GTP encounter complexes, the MSM accurately discovers the native pose of GTP at its designated catalytic site in crystallographic precision. However, the series of events exhibit signatures of conformational plasticity in which the protein remains trapped in multiple non-native conformations even when GTP has already located itself in its native binding site. The investigation demonstrates mechanistic relays pertaining to simultaneous fluctuations of switch 1 and switch 2 residues which remain most instrumental in maneuvering the GTP-binding process. Scanning of the crystallographic database reveals close resemblance between observed non-native GTP binding poses and precedent crystal structures of substrate-bound GTPase, suggesting potential roles of these binding-competent intermediates in allosteric regulation of the recognition process.

A deep encoder-decoder framework for identifying distinct ligand binding pathways

The pathway(s) that a ligand would adopt en route to its trajectory to the native pocket of the receptor protein act as a key determinant of its biological activity. While Molecular Dynamics (MD) simulations have emerged as the method of choice for modeling protein-ligand binding events, the high dimensional nature of the MD-derived trajectories often remains a barrier in the statistical elucidation of distinct ligand binding pathways due to the stochasticity inherent in the ligand's fluctuation in the solution and around the receptor. Here, we demonstrate that an autoencoder based deep neural network, trained using an objective input feature of a large matrix of residue-ligand distances, can efficiently produce an optimal low-dimensional latent space that stores necessary information on the ligand-binding event. In particular, for a system of L99A mutant of T4 lysozyme interacting with its native ligand, benzene, this deep encoder-decoder framework automatically identifies multiple distinct recognition pathways, without requiring user intervention. The intermediates involve the spatially discrete location of the ligand in different helices of the protein before its eventual recognition of native pose. The compressed subspace derived from the autoencoder provides a quantitatively accurate measure of the free energy and kinetics of ligand binding to the native pocket. The investigation also recommends that while a linear dimensional reduction technique, such as time-structured independent component analysis, can do a decent job of state-space decomposition in cases where the intermediates are long-lived, autoencoder is the method of choice in systems where transient, low-populated intermediates can lead to multiple ligand-binding pathways.

Resolving Protein Conformational Plasticity and Substrate Binding via Machine Learning

A long-standing target in elucidating the biomolecular recognition process is the identification of binding-competent conformations of the receptor protein. However, protein conformational plasticity and the stochastic nature of the recognition processes often preclude the assignment of a specific protein conformation to an individual ligand-bound pose. Here, we demonstrate that a computational framework coined as RF-TICA-MD, which integrates an ensemble decision-tree-based Random Forest (RF) machine learning (ML) technique with an unsupervised dimension reduction approach time-structured independent component analysis (TICA), provides an efficient and unambiguous solution toward resolving protein conformational plasticity and the substrate binding process. In particular, we consider multimicrosecond-long molecular dynamics (MD) simulation trajectories of a ligand recognition process in solvent-inaccessible cavities of archetypal proteins T4 lysozyme and cytochrome P450cam. We show that in a scenario in which clear correspondence between protein conformation and binding-competent macrostates could not be obtained via an unsupervised dimension reduction approach, an a priori decision-tree-based supervised classification of the simulated recognition trajectories via RF would help characterize key amino acid residue pairs of the protein that are deemed sensitive for ligand binding. A subsequent unsupervised dimensional reduction of the selected residue pairs via TICA would then delineate a conformational landscape of protein which is able to demarcate ligand-bound poses from unbound ones. The proposed RF-TICA-MD approach is shown to be data agnostic and found to be robust when using other ML-based classification methods such as XGBoost. As a promising spinoff of the protocol, the framework is found to be capable of identifying distal protein locations which would be allosterically important for ligand binding and would characterize their roles in recognition pathways. A Python implementation of a proposed ML workflow is available in GitHub https://github.com/navjeet0211/rf-tica-md.

Long-time-step molecular dynamics can retard simulation of protein-ligand recognition process

Molecular dynamics (MD) simulation of biologically relevant processes at realistic time scale and atomistic precision is generally limited by prohibitively large computational cost, due to its restriction of using an ultrashort integration time step (1-2 fs). A popular numerical recipe to reduce the associated computational burden is adopting schemes that would allow relatively longer-time-step for MD propagation. Here, we explore the perceived potential of one of the most frequently used long-time-step protocols, namely the hydrogen mass repartitioning (HMR) approach, in alleviating the computational overhead associated with simulation of the kinetic process of protein-ligand recognition events. By repartitioning the mass of heavier atoms to their linked hydrogen atoms, HMR leverages around twofold longer time step than regular simulation, holding promise of significant performance boost. However, our probe into direct simulation of the protein-ligand recognition event, one of the computationally most challenging processes, shows that long-time-step HMR MD simulations do not necessarily translate to a computationally affordable solution. Our investigations spanning cumulative 176 μs in three independent proteins (T4 lysozyme, sensor domain of MopR, and galectin-3) show that long-time-step HMR-based MD simulations can catch the ligand in its act of recognizing the native cavity. But, as a major caveat, the ligand is found to require significantly longer time to identify buried native protein cavity in an HMR MD simulation than regular simulation, thereby defeating the purpose of its usage for performance upgrade. A molecular analysis shows that the longer time required by a ligand to recognize the protein in HMR is rooted in faster diffusion of the ligand, which reduces the survival probability of decisive on-pathway metastable intermediates, thereby slowing down the eventual recognition process at the native cavity. Together, the investigation stresses careful assessment of pitfalls of long-time-step algorithms before attempting to utilize them for higher performance for biomolecular recognition simulations.

Markov State Models Underlying the N-Terminal Premodel of TOPK/PBK

Lymphokine-activated killer T-cell-originated protein kinase (TOPK) is a potential target for cancer therapy. To explore the micromechanism, we proposed the N-terminal premodel (NTPM) of the TOPK monomer via homology modeling and molecular dynamic simulations and analyzed the conformational dynamics by Markov state model analysis. The electronegative insert (ENI) motif of the NTPM can be opened with a small probability under wild type, regulated by the so-called "N-C" interaction zone consisting of the N-terminal head, the coil between β₃-strand and αC-helix, and the ENI motif. Glutamate substitution at threonine residue 9 or tyrosine residue 74 promotes the closed-open transition, revealing the details of phosphorylation. Allosteric effects induce functionally relevant structural changes, such as increased structural flexibility and active sites, which are thought to be necessary for further activation or binding. These findings provide rational structural templates for designing state-dependent inhibitors and give insight into the molecular regulatory mechanisms of TOPK monomers.

Solvent's Role in Cavity-Ligand Recognition Would Depend on the Mode of Ligand Diffusion

Solvent is known to play crucial roles in dictating the thermodynamics and kinetics of the biomolecular recognition process. Here, we show that the extent of significance of water in modulating the ligand recognition process is critically contingent on the ligand diffusion and on the constraints introduced on it. Toward the end, we use a well-known prototypical system of spherical ligand diffusing freely toward a hydrophobic concave cavity in explicit water. We analyze a large series of adaptively sampled unbiased molecular dynamics simulation trajectories within the framework of time-structured independent component analysis (TICA). Our quantitative investigations reveal that water would play a significant role in the ligand recognition process, provided that the ligand is constricted to diffuse along a centro-symmetric fashion. On the contrary, water's contribution in the ligand recognition process would diminish to a negligible value if the ligand freely diffuses toward the pocket. A Markov state model (MSM) constructed using the simulated trajectories identifies a set of transiently populated metastable states comprising partially ligand-unbound macro states, alongside ligand-bound and ligand-unbound pose and gives rise to multiple transition paths of ligand in its way to the hydrophobic cavity. Lifting the restriction on ligand movement changes its binding pathway, time scales, and the extent of the role of solvent in modulating the recognition process.

Molecular dynamics simulations elucidate oligosaccharide recognition pathways by galectin-3 at atomic resolution

The recognition of carbohydrates by lectins plays key roles in diverse cellular processes such as cellular adhesion, proliferation, and apoptosis, which makes it a therapeutic target of significance against cancers. One of the most functionally active lectins, galectin-3 is distinctively known for its specific binding affinity toward β-galactoside. However, despite the prevalence of high-resolution crystallographic structures, the mechanistic basis and more significantly, the dynamic process underlying carbohydrate recognition by galectin-3 are currently elusive. To this end, we employed extensive Molecular Dynamics simulations to unravel the complete binding event of human galectin-3 with its native natural ligand N-acetyllactosamine (LacNAc) at atomic precision. The simulation trajectory demonstrates that the oligosaccharide diffuses around the protein and eventually identifies and binds to the biologically designated binding site of galectin-3 in real time. The simulated bound pose correlates with the crystallographic pose with atomic-level accuracy and recapitulates the signature stabilizing galectin-3/oligosaccharide interactions. The recognition pathway also reveals a set of transient non-native ligand poses in its course to the receptor. Interestingly, kinetic analysis in combination with a residue-level picture revealed that the key to the efficacy of a more active structural variant of the LacNAc lay in the ligand's resilience against disassociation from galectin-3. By catching the ligand in the act of finding its target, our investigations elucidate the detailed recognition mechanism of the carbohydrate-binding domain of galectin-3 and underscore the importance of ligand-target binary complex residence time in understanding the structure-activity relationship of cognate ligands.

Molecular Mechanism of Methanol Inhibition in CALB-Catalyzed Alcoholysis: Analyzing Molecular Dynamics Simulations by a Markov State Model

Lipases are widely used enzymes that catalyze hydrolysis and alcoholysis of fatty acid esters. At high concentrations of small alcohols such as methanol or ethanol, many lipases are inhibited by the substrate. The molecular basis of the inhibition of Candida antarctica lipase B (CALB) by methanol was investigated by unbiased molecular dynamics (MD) simulations, and the substrate binding kinetics was analyzed by Markov state models (MSMs). The modeled fluxes of productive methanol binding at concentrations between 50 mM and 5.5 M were in good agreement with the experimental activity profile of CALB, with a peak at 300 mM. The kinetic and structural analysis uncovered the molecular basis of CALB inhibition. Beyond 300 mM, the kinetic bottleneck results from crowding of methanol in the substrate access channel, which is caused by the gradual formation of methanol patches close to Leu140 (helix α5), Leu278, and Ile285 (helix α10) at a distance of 4-5 Å from the active site. Our findings demonstrate the usefulness of unbiased MD simulations to study enzyme-substrate interactions at realistic substrate concentrations and the feasibility of scale-bridging by an MSM analysis to derive kinetic information.

A deep autoencoder framework for discovery of metastable ensembles in biomacromolecules

Biomacromolecules manifest dynamic conformational fluctuation and involve mutual interconversion among metastable states. A robust mapping of their conformational landscape often requires the low-dimensional projection of the conformational ensemble along optimized collective variables (CVs). However, the traditional choice for the CV is often limited by user-intuition and prior knowledge about the system, and this lacks a rigorous assessment of their optimality over other candidate CVs. To address this issue, we propose an approach in which we first choose the possible combinations of inter-residue Cα-distances within a given macromolecule as a set of input CVs. Subsequently, we derive a non-linear combination of latent space embedded CVs via auto-encoding the unbiased molecular dynamics simulation trajectories within the framework of the feed-forward neural network. We demonstrate the ability of the derived latent space variables in elucidating the conformational landscape in four hierarchically complex systems. The latent space CVs identify key metastable states of a bead-in-a-spring polymer. The combination of the adopted dimensional reduction technique with a Markov state model, built on the derived latent space, reveals multiple spatially and kinetically well-resolved metastable conformations for GB1 β-hairpin. A quantitative comparison based on the variational approach-based scoring of the auto-encoder-derived latent space CVs with the ones obtained via independent component analysis (principal component analysis or time-structured independent component analysis) confirms the optimality of the former. As a practical application, the auto-encoder-derived CVs were found to predict the reinforced folding of a Trp-cage mini-protein in aqueous osmolyte solution. Finally, the protocol was able to decipher the conformational heterogeneities involved in a complex metalloenzyme, namely, cytochrome P450.

Development of MEMS directed evolution strategy for multiplied throughput and convergent evolution of cytochrome P450 enzymes

Directed evolution (DE) inspired by natural evolution (NE) has been achieving tremendous successes in protein/enzyme engineering. However, the conventional "one-protein-for-one-task" DE cannot match the "multi-proteins-for-multi-tasks" NE in terms of screening throughput and efficiency, thus often failing to meet the fast-growing demands for biocatalysts with desired properties. In this study, we design a novel "multi-enzymes-for-multi-substrates" (MEMS) DE model and establish the proof-of-concept by running a NE-mimicking and higher-throughput screening on the basis of "two-P450s-against-seven-substrates" (2P×7S) in one pot. With the multiplied throughput and improved hit rate, we witness a series of convergent evolution events of the two archetypal cytochrome P450 enzymes (P450 BM3 and P450cam) in laboratory. It is anticipated that the new strategy of MEMS DE will find broader application for a larger repertoire of enzymes in the future. Furthermore, structural and substrate docking analysis of the two functionally convergent P450 variants provide important insights into how distinct P450 active-sites can reach a common catalytic goal.

Reconciling conformational heterogeneity and substrate recognition in cytochrome P450

Cytochrome P450, the ubiquitous metalloenzyme involved in detoxification of foreign components, has remained one of the most popular systems for substrate-recognition process. However, despite being known for its high substrate specificity, the mechanistic basis of substrate-binding by archetypal system cytochrome P450cam has remained at odds with the contrasting reports of multiple diverse crystallographic structures of its substrate-free form. Here, we address this issue by elucidating the probability of mutual dynamical transition to the other crystallographic pose of cytochrome P450cam and vice versa via unbiased all-atom computer simulation. A robust Markov state model, constructed using adaptively sampled 84-μs-long molecular dynamics simulation trajectories, maps the broad and heterogenous P450cam conformational landscape into five key substates. In particular, the Markov state model identifies an intermediate-assisted dynamic equilibrium between a pair of conformations of P450cam, in which the substrate-recognition sites remain "closed" and "open," respectively. However, the estimate of a significantly higher stationary population of closed conformation, coupled with faster rate of open → closed transition than its reverse process, dictates that the net conformational equilibrium would be swayed in favor of "closed" conformation. Together, the investigation quantitatively infers that although a potential substrate of cytochrome P450cam would, in principle, explore a diverse array of conformations of substrate-free protein, it would mostly encounter a "closed" or solvent-occluded conformation and hence would follow an induced-fit-based recognition process. Overall, the work reconciles multiple precedent crystallographic, spectroscopic investigations and establishes how a statistical elucidation of conformational heterogeneity in protein would provide crucial insights in the mechanism of potential substrate-recognition process.

Conformational plasticity and dynamic interactions of the N-terminal domain of the chemokine receptor CXCR1

The dynamic interactions between G protein-coupled receptors (GPCRs) and their cognate protein partners are central to several cell signaling pathways. For example, the association of CXC chemokine receptor 1 (CXCR1) with its cognate chemokine, interleukin-8 (IL8 or CXCL8) initiates pathways leading to neutrophil-mediated immune responses. The N-terminal domain of chemokine receptors confers ligand selectivity, but unfortunately the conformational dynamics of this intrinsically disordered region remains unresolved. In this work, we have explored the interaction of CXCR1 with IL8 by microsecond time scale coarse-grain simulations, complemented by atomistic models and NMR chemical shift predictions. We show that the conformational plasticity of the apo-receptor N-terminal domain is restricted upon ligand binding, driving it to an open C-shaped conformation. Importantly, we corroborated the dynamic complex sampled in our simulations against chemical shift perturbations reported by previous NMR studies and show that the trends are similar. Our results indicate that chemical shift perturbation is often not a reporter of residue contacts in such dynamic associations. We believe our results represent a step forward in devising a strategy to understand intrinsically disordered regions in GPCRs and how they acquire functionally important conformational ensembles in dynamic protein-protein interfaces.

Pervasive cooperative mutational effects on multiple catalytic enzyme traits emerge via long-range conformational dynamics

Multidimensional fitness landscapes provide insights into the molecular basis of laboratory and natural evolution. To date, such efforts usually focus on limited protein families and a single enzyme trait, with little concern about the relationship between protein epistasis and conformational dynamics. Here, we report a multiparametric fitness landscape for a cytochrome P450 monooxygenase that was engineered for the regio- and stereoselective hydroxylation of a steroid. We develop a computational program to automatically quantify non-additive effects among all possible mutational pathways, finding pervasive cooperative signs and magnitude epistasis on multiple catalytic traits. By using quantum mechanics and molecular dynamics simulations, we show that these effects are modulated by long-range interactions in loops, helices and β-strands that gate the substrate access channel allowing for optimal catalysis. Our work highlights the importance of conformational dynamics on epistasis in an enzyme involved in secondary metabolism and offers insights for engineering P450s.

Connecting conformational dynamics and epistasis has so far been limited to a few proteins and a single fitness trait. Here, the authors provide evidence of positive epistasis on multiple catalytic traits in the evolution and dynamics of engineered cytochrome P450 monooxygenase, offering insights for in silico protein design.

An Appraisal of Computer Simulation Approaches in Elucidating Biomolecular Recognition Pathways

Computer simulation approaches in biomolecular recognition processes have come a long way. In this Perspective, we highlight a series of recent success stories in which computer simulations have played a remarkable role in elucidating the atomic resolution mechanism of kinetic processes of protein-ligand binding in a quantitative fashion. In particular, we show that a robust combination of unbiased simulation, harnessed by a high-fidelity computing environment, and Markov state modeling approaches has been instrumental in revealing novel protein-ligand recognition pathways in multiple systems. We also elucidate the role of recent developments in enhanced sampling approaches in providing the much-needed impetus in accelerating simulation of the ligand recognition process. We identify multiple key issues, including force fields and the sampling bottleneck, which are currently preventing the field from achieving quantitative reconstruction of experimental measurements. Finally, we suggest a possible way forward via adoption of multiscale approaches and coarse-grained simulations as next steps toward efficient elucidation of ligand binding kinetics.

Simulating Multiple Substrate-Binding Events by γ-Glutamyltransferase Using Accelerated Molecular Dynamics

γ-Glutamyltransferase (GGT) is an enzyme that uses γ-glutamyl compounds as substrates and catalyzes their transfer to a water molecule or an acceptor substrate with varied physiological function in bacteria, plants, and animals. Crystal structures of GGT are known for different species and in different states of the chemical reaction; however, the structural dynamics of the substrate binding to the catalytic site of GGT are unknown. Here, we modeled Escherichia coli GGT's glutamine binding by using a swarm of accelerated molecular dynamics (aMD) simulations. Characterization of multiple binding events identified three structural binding motifs composed of polar residues in the binding pocket that govern glutamine binding into the active site. Simulated open and closed conformations of a lid-loop protecting the binding cavity suggest its role as a gating element by allowing or blocking substrates entry into the binding pocket. Partially open states of the lid-loop are accessible within thermal fluctuations, while the estimated free energy cost of a complete open state is 2.4 kcal/mol. Our results suggest that both specific electrostatic interactions and GGT conformational dynamics dictate the molecular recognition of substrate-GGT complexes.

Capturing Protein-Ligand Recognition Pathways in Coarse-Grained Simulation

Protein-ligand recognition is dynamic and complex. A key approach in deciphering the mechanism underlying the recognition process is to capture the kinetic process of the ligand in its act of binding to its designated protein cavity. Toward this end, ultralong all-atom molecular dynamics simulation has recently emerged as a popular method of choice because of its ability to record these events at high spatial and temporal resolution. However, success via this route comes at an exorbitant computational cost. Herein, we demonstrate that coarse-grained models of the protein, when systematically optimized to maintain its tertiary fold, can capture the complete process of spontaneous protein-ligand binding from bulk media to the cavity at crystallographic precision and within wall clock time that is orders of magnitude shorter than that of all-atom simulations. The exhaustive sampling of ligand exploration in protein and solvent, harnessed by coarse-grained simulation, leads to elucidation of new ligand recognition pathways and discovery of non-native binding poses.

Global and Kinetic Profiles of Substrate Diffusion in Candida antarctica Lipase B: Molecular Dynamics with the Markov-State Model

Profiling substrate diffusion pathways with kinetic information, which accounts for the dynamic nature of enzyme-substrate interaction, can enable molecular reengineering of enzymes and process optimization of enzymatic catalysis. Candida antarctica lipase B (CALB) is extensively used for producing various chemicals because of its rich catalytic mechanisms, broad substrate spectrum, thermal stability, and tolerance to organic solvents. In this study, an all-atom molecular dynamics (MD) combined with Markov-state models (MSMs) implemented in pyEMMA was proposed to simulate diffusion pathways of 4-nitrophenyl ester (4NPE), a commonly used substrate, from the surface into the active site of CALB. Six important metastable conformations of CALB were identified in the diffusion process, including a closed state. An induced-fit mechanism incorporating multiple pathways with molecular information was proposed, which might find unprecedented applications for the rational design of lipase for green catalysis.

Spontaneous Ligand Access Events to Membrane-Bound Cytochrome P450 2D6 Sampled at Atomic Resolution

The membrane-anchored enzyme Cytochrome P450 2D6 (CYP2D6) is involved in the metabolism of around 25% of marketed drugs and its metabolic performance shows a high interindividual variation. While it was suggested that ligands access the buried active site of the enzyme from the membrane, no proof from unbiased simulations has been provided to support this hypothesis. Laboratory experiments fail to capture the access process which is suspected to influence binding kinetics. Here, we applied unbiased molecular dynamics (MD) simulations to investigate the access of ligands to wild-type CYP2D6, as well as the allelic variant CYP2D6*53. In multiple simulations, substrates accessed the active site of the enzyme from the protein-membrane interface to ultimately adopt a conformation that would allow a metabolic reaction. We propose the necessary steps for ligand access and the results suggest that the increased metabolic activity of CYP2D6*53 might be caused by a facilitated ligand uptake.

In Silico Reoptimization of Binding Affinity and Drug-Resistance Circumvention Ability in Kinase Inhibitors: A Case Study with RL-45 and Src Kinase

A major bottleneck in the development of kinase inhibitors has been the onset of drug resistance around the gatekeeper residues of Src kinase. Although recent times have seen the reports of certain second-generation kinase inhibitors which are capable of bypassing the drug resistance by circumventing kinase mutation, their kinase-binding efficacy has remained considerably weaker than that of the classical adenosine 5'-triphosphate-competitive kinase inhibitors. Using a recently synthesized second-generation kinase inhibitor RL-45 as a template, the current work integrates fragment-based drug discovery and quantitative structure-activity relationship study with enhanced molecular dynamics simulation approaches, namely, metadynamics and replica exchange free-energy perturbation, and demonstrates how one can optimally redesign and assess novel Src kinase inhibitors, by minimal introduction of new functional moieties around template kinase inhibitor. Interestingly, unlike many synthetic kinase inhibitors, these in silico optimized small-molecule derivatives of RL-45 are found to be potentially capable of serving dual purposes, crucial for efficacy of an ideal kinase inhibitor: (a) circumventing gatekeeper residue mutation-related drug resistance in Src kinase, unlike many commercial kinase inhibitors and (b) manifesting superior resilience against unbinding from the kinase active site. The computer simulation, boosted by enhanced sampling techniques, further reveals that these designed inhibitors bring about key interactions in the form of significantly long-standing hydrogen bonds and hydrophobic pocket otherwise weak in the template bioactive kinase inhibitor, which enhance the binding efficacy of these newly designed ligands in the kinase-binding pocket.