Unbinding Process Research Articles

Decoding Inhibitor Egression from Wild-Type and G2019S Mutant LRRK2 Kinase: Insights into Unbinding Mechanisms for Precision Drug Design in Parkinson's Disease.

Leucine-rich repeat kinase 2 (LRRK2) remains a viable target for drug development since the discovery of the association of its mutations with Parkinson's disease (PD). G2019S (in the kinase domain) is the most common mutation for LRRK2-based PD. Though various types of inhibitors have been developed for the kinase domain to reduce the effect of the mutation, understanding the working of these inhibitors at the molecular level is still ongoing. This study focused on the exploration of the dissociation mechanism (pathways) of inhibitors from (WT and G2019S) LRRK2 kinase (using homology model CHK1 kinase), which is one of the crucial aspects in drug discovery. Here, two ATP-competitive type I inhibitors, PF-06447475 and MLi-2 (Comp1 and Comp2 ), and one non-ATP-competitive type II inhibitor, rebastinib (Comp3), were considered for this investigation. To study the unbinding process, random accelerated molecular dynamics simulations were performed. The binding free energies of the three inhibitors for different egression paths were determined using umbrella sampling. This work found four major egression pathways that were adopted by the inhibitors Comp1 (path1, path2, and path3), Comp2 (path1, path2 and path3), and Comp3 (path3 and path4). Also, the mechanism of unbinding for each path and key residues involved in unbinding were explored. Mutation was not observed to impact the preference of the particular egression pathways for both LRRK2-Comp1 and -Comp2 systems. However, the findings suggested that the size of the inhibitor molecules might have an effect on the preference of the egression pathways. The binding energy and residence time of the inhibitors followed a similar trend to experimental observations. The findings of this work might provide insight into designing more potent inhibitors for the G2019S LRRK2 kinase.

Role of UDP-N-acetylmuramic acid in the regulation of MurA activity revealed by molecular dynamics simulations.

The peptidoglycan biosynthesis pathway plays a vital role in bacterial cells, and facilitates peptidoglycan layer formation, a fundamental structural component of the bacterial cell wall. The enzymes in this pathway are candidates for antibiotic development, as most do not have mammalian homologues. The UDP-N-acetylglucosamine (UNAG) enolpyruvyl transferase enzyme (MurA) in the peptidoglycan pathway cytoplasmic step is responsible for the phosphoenolpyruvate (PEP)-UNAG catalytic reaction, forming UNAG enolpyruvate and inorganic phosphate. Reportedly, UDP-N-acetylmuramic acid (UNAM) binds tightly to MurA forming a dormant UNAM-PEP-MurA complex and acting as a MurA feedback inhibitor. MurA inhibitors are complex, owing to competitive binding interactions with PEP, UNAM, and UNAG at the MurA active site. We used computational methods to explore UNAM and UNAG binding. UNAM showed stronger hydrogen-bond interactions with the Arg120 and Arg91 residues, which help to stabilize the closed conformation of MurA, than UNAG. Binding free energy calculations using end-point computational methods showed that UNAM has a higher binding affinity than UNAG, when PEP is attached to Cys115. The unbinding process, simulated using τ-random acceleration molecular dynamics, showed that UNAM has a longer relative residence time than UNAG, which is related to several complex dissociation pathways, each with multiple intermediate metastable states. This prevents the loop from opening and exposing the Arg120 residue to accommodate UNAG and potential new ligands. Moreover, we demonstrate the importance of Cys115-linked PEP in closed-state loop stabilization. We provide a basis for evaluating novel UNAM analogues as potential MurA inhibitors. PUBLIC SIGNIFICANCE: MurA is a critical enzyme involved in bacterial cell wall biosynthesis and is involved in antibiotic resistance development. UNAM can remain in the target protein's active site for an extended time compared to its natural substrate, UNAG. The prolonged interaction of this highly stable complex known as the 'dormant complex' comprises UNAM-PEP-MurA and offers insights into antibiotic development, providing potential options against drug-resistant bacteria and advancing our understanding of microbial biology.

Exploring Product Release from Yeast Cytosine Deaminase with Metadynamics.

The yeast cytosine deaminase (yCD) enzyme/5-fluorocytosine prodrug system is a promising candidate for targeted chemotherapeutics. After conversion of the prodrug into the toxic chemotherapeutic drug, 5-fluorouracil (5-FU), the slow product release from the enzyme limits the overall catalytic efficiency of the enzyme/prodrug system. Here, we present a computational study of the product release of the anticancer drug, 5-FU, from yCD using metadynamics. We present a comparison of the 5-FU drug to the natural enzyme product, uracil. We use volume-based metadynamics to compute the free energy landscape for product release and show a modest binding affinity for the product to the enzyme, consistent with experiments. Next, we use infrequent metadynamics to estimate the unbiased release rate from Kramers time-dependent rate theory and find a favorable comparison to experiment with a slower rate of product release for the 5-FU system. Our work demonstrates how adaptive sampling methods can be used to study the protein-ligand unbinding process for engineering enzyme/prodrug systems and gives insights into the molecular mechanism of product release for the yCD/5-FU system.

The mutual and dynamic role of TSPO and ligands in their binding process: An example with PK-11195

Translocator protein (TSPO) is an 18 kDa transmembrane protein, localized primarily on the outer mitochondrial membrane. It has been found to be involved in various physiological processes and pathophysiological conditions. Though studies on its structure have been performed only recently, there is little information on the nature of dynamics and doubts about some structures referenced in the literature, especially the NMR structure of mouse TSPO. In the present work, we thoroughly study the dynamics of mouse TSPO protein by means of atomistic molecular dynamics simulations, in presence as well as in absence of the diagnostic ligand PKA. We considered two starting structures: the NMR structure and a homology model (HM) generated on the basis of X-ray structures from bacterial TSPO. We examine the conformational landscape in both the modes for both starting points, in presence and absence of the ligand, in order to measure its impact for both structures. The analysis highlights high flexibility of the protein globally, but NMR simulations show a surprisingly flexibility even in the presence of the ligand. Interestingly, this is not the case for HM calculations, to the point that the ligand seems not so stable as in the NMR system and an unbinding event process is partially sampled. All those results tend to show that the NMR structure of mTSPO seems not deficient but is just in another portion of the global conformation space of TSPO.

Exploring transmembrane allostery in the MexB: DB08385 variant as a promising inhibitor-like candidate against Pseudomonas aeruginosa antibiotic resistance: a computational study.

Pseudomonas aeruginosa, a formidable pathogen renowned for its antimicrobial resistance, poses a significant threat to immunocompromised individuals. In this regard, the MexAB-OprM efflux pump acts as a pivotal line of defense by extruding antimicrobials from bacterial cells. The inner membrane homotrimeric protein MexB captures antibiotics and translocates them into the outer membrane OprM channel protein connected through the MexA adaptor protein. Despite extensive efforts, competitive inhibitors targeting the tight (T) protomer of the MexB protein have not received FDA approval for medical use. Over the past few years, allosteric inhibitors have become popular as alternatives to the classical competitive inhibitor-based approach because of their higher specificity, lower dosage, and reduced toxicological effects. Hence, in this study, we unveiled the existence of a transmembrane allosteric binding pocket of MexB inspired by the recent discovery of an important allosteric inhibitor, BDM88855, for the homolog AcrB protein. While repurposing BDM88855 proved ineffective in controlling the MexB loose (L) protomer, our investigation identified a promising alternative: a chlorine-containing variant of DB08385 (2-Cl DB08385 or Variant 1). Molecular dynamics simulations, including binding free energy estimation coupled with heterogeneous dielectric implicit membrane model (implicit-membrane MM/PBSA), interaction entropy (IE) analysis and potential of mean force (PMF) calculation, demonstrated Variant 1's superior binding affinity to the transmembrane pocket, displaying the highest energy barrier in the ligand unbinding process. To elucidate the allosteric crosstalk between the transmembrane and porter domain of MexB, we employed the 'eigenvector centrality' measure in the linear mutual information obtained from the protein correlation network. Notably, this study confirmed the presence of an allosteric transmembrane site in the MexB L protomer. In addition to this, Variant 1 emerged as a potent regulator of allosteric crosstalk, inducing an 'O-L intermediate state' in the MexB L protomer. This induced state might hold the potential to diminish substrate intake into the access pocket, leading to the ineffective efflux of antibiotics.

Rate limiting step of the allosteric activation of the bacterial adhesin FimH investigated by molecular dynamics simulations.

The bacterial adhesin FimH is a model for the study of protein allostery because its structure has been resolved in multiple configurations, including the active and the inactive state. FimH consists of a pilin domain (PD) that anchors it to the rest of the fimbria and an allosterically regulated lectin domain (LD) that binds mannose on the surface of infected cells. Under normal conditions, the two domains are docked to each other and LD binds mannose weakly. However, in the presence of tensile force generated by shear the domains separate and conformational changes propagate across LD resulting in a stronger bond to mannose. Recently, the crystallographic structure of a variant of FimH has been resolved, called , where PD contains 10 mutations near the inter-domain interface. Although the X-ray structures of FimH and are almost identical, experimental evidence shows that is activated even in the absence of shear. Here, molecular dynamics simulations combined with the Jarzynski equality were used to investigate the discrepancy between the crystallographic structures and the functional assays. The results indicate that the free energy barrier of the unbinding process between LD and PD is drastically reduced in . Rupture of inter-domain hydrogen bonds involving R166 constitutes a rate limiting step of the domain separation process and occurs more readily in than FimH. In conclusion, the mutations in shift the equilibrium toward an equal occupancy of bound and unbound states for LD and PD by reducing a rate limiting step.

Structural insights into the lead identification of sub-type selective PDE4B inhibitors from plant bioactive molecule analogues

Phosphodiesterase-4B (PDE4B), a sub-type of PDE4, is highly expressed in the brain and associated with schizophrenia, neuroinflammation, and cognition. A selective inhibitor of PDE4B is anticipated to be a promising drug candidate for addressing related disorders. It is crucial to identify the factors that influence inhibitor specificity and the binding process to design PDE4B inhibitors with high efficiency and specificity. This includes investigating the role of various essential residues within the active site pockets (Q, M, and S). Therefore, molecular docking and atomistic molecular dynamics simulations were performed to assess the inhibitory potential of benzosuberene-sulfone (BSS) molecules against PDE4B compared to the standard drug roflumilast. The S6 and S10 BSS analogues formed stable complexes and exhibited lower binding free energy than the standard molecule. Additionally, it was observed that the S6 and S10 accommodated well inside the Q and S-pockets of PDE4B. The study also highlighted the importance of amino acid residue Met-431 (Q-pocket residue) to achieve sub-type selectivity. Furthermore, the pulling simulations revealed that the S6 molecule required the highest force to unbind from the binding cavity and showed key interactions with Met-431 during the unbinding process. The umbrella sampling simulations revealed binding free energies of -69.09 ± 1.37 kJ/mol for S6 and -55.11 ± 1.69 kJ/mol for roflumilast, underscoring a higher affinity of S6 for PDE4B. This study presented compelling evidences that the S6 molecule could be a promising lead for the selective inhibition of PDE4B.

Lessons learned through weighted ensemble simulations of a (surprisingly) challenging protein-ligand unbinding system.

A major challenge in applying the weighted ensemble enhanced sampling strategy to simulate a rare-event process of interest is to identify an effective progress coordinate. Furthermore, the algorithm for the resampling procedure of the strategy can impact both the number and diversity of representative pathways generated. Here, I will share lessons learned from the application of atomistic, weighted ensemble simulations to a surprisingly challenging protein-ligand unbinding process that features a highly-charged ligand. These lessons include key considerations regarding the choice of progress coordinate and resampling algorithm, including binless resamplers.

Unraveling the Abnormal Molecular Mechanism of Suicide Inhibition of Cytochrome P450 3A4.

Suicide inhibition of the CYP3A4 enzyme by a drug inactivates the enzyme in the drug biotransformation process and often shows safety concerns about the drug. Despite extensive experimental studies, the abnormal molecular mechanism of a suicide inhibitor that forms a covalent bond with the residue far away from the catalytically active center of CYP3A4 inactivating the enzyme remains elusive. Here, the authors used molecular simulation approaches to study in detail how diquinone methide (DQR), the metabolite product of raloxifene, unbinds from CYP3A4 and inactivates the enzyme at the atomistic level. The results clearly indicate that in one of the intermediate states formed in its unbinding process, DQR covalently binds to Cys239, a residue far away from the catalytically active center of CYP3A4, and hinders the substrate from entering or leaving the enzyme. This work therefore provides an unprecedented way of clarifying the abnormal mechanism of suicide inhibition of the CYP3A4 enzyme.

Water Network in the Binding Pocket of Fluorinated BPTI-Trypsin Complexes─Insights from Simulation and Experiment.

Structural waters in the S1 binding pocket of β-trypsin are critical for the stabilization of the complex of β-trypsin with its inhibitor bovine pancreatic trypsin inhibitor (BPTI). The inhibitor strength of BPTI can be modulated by replacing the critical lysine residue at the P1 position by non-natural amino acids. We study BPTI variants in which the critical Lys15 in BPTI has been replaced by α-aminobutyric acid (Abu) and its fluorinated derivatives monofluoroethylglycine (MfeGly), difluoroethylglycine (DfeGly), and trifluoroethylglycine (TfeGly). We investigate the hypothesis that additional water molecules in the binding pocket can form specific noncovalent interactions with the fluorinated side chains and thereby act as an extension of the inhibitors. We report potentials of mean force (PMF) of the unbinding process for all four complexes and enzyme activity inhibition assays. Additionally, we report the protein crystal structure of the Lys15MfeGly-BPTI-β-trypsin complex (pdb: 7PH1). Both experimental and computational data show a stepwise increase in inhibitor strength with increasing fluorination of the Abu side chain. The PMF additionally shows a minimum for the encounter complex and an intermediate state just before the bound state. In the bound state, the computational analysis of the structure and dynamics of the water molecules in the S1 pocket shows a highly dynamic network of water molecules that does not indicate a rigidification or stabilizing trend in regard to energetic properties that could explain the increase in inhibitor strength. The analysis of the energy and the entropy of the water molecules in the S1 binding pocket using grid inhomogeneous solvation theory confirms this result. Overall, fluorination systematically changes the binding affinity, but the effect cannot be explained by a persistent water network in the binding pocket. Other effects, such as the hydrophobicity of fluorinated amino acids and the stability of the encounter complex as well as the additional minimum in the potential of mean force in the bound state, likely influence the affinity more directly.

Comparison of the unbinding process of RBD-ACE2 complex between SARS-CoV-2 variants (Delta, delta plus, and Lambda): A steered molecular dynamics simulation

ABSTRACT The Delta and Lambda SARS-CoV-2 variants are characterised with the double mutations T478K/L452R, and F490S/L452Q, respectively, in their receptor binding domains (RBDs) of the spike proteins. The Delta Plus variant is one of the Delta variant types with K417N substitution in RBD. Basically, SARS-CoV-2 variants need to bind to Angiotensin-Converting Enzyme 2(ACE2) to enter the host cell. Here, the effects of these mutations were studied on the binding affinity of RBD-ACE2 using steered molecular dynamics (SMD) simulation. According to the results, the rupture force of RBD-ACE2 for the wild type and Delta variant were close to each other, and the minimum and maximum rupture forces occurred in the Delta Plus and Lambda variants, respectively. While T478K and L452R occur adjacent to the RBD binding site and have no major effect on the rupture force of the RBD-ACE2 complex, K417N occurs exactly at the binding site of RBD-ACE2 and causes a decrease in the rupture force for the Delta Plus variant. The results of this study showed that T478K, L452R, and K417N could not justify the infection and transmissibility of the Delta and Delta Plus variants. However, L452Q and F490S were able to affect the binding affinity of the RBD-ACE2 complex.

Comparative study of the unbinding process of some HTLV-1 protease inhibitors using unbiased molecular dynamics simulations.

The HTLV-1 protease is one of the major antiviral targets to overwhelm this virus. Several research groups have developed protease inhibitors, but none has been successful. In this regard, developing new HTLV-1 protease inhibitors to fix the defects in previous inhibitors may overcome the lack of curative treatment for this oncovirus. Thus, we decided to study the unbinding pathways of the most potent (compound 10, PDB ID 4YDF, Ki = 15 nM) and one of the weakest (compound 9, PDB ID 4YDG, Ki = 7900 nM) protease inhibitors, which are very structurally similar. We conducted 12 successful short and long simulations (totaling 14.8 μs) to unbind the compounds from two monoprotonated (mp) forms of protease using the Supervised Molecular Dynamics (SuMD) without applying any biasing force. The results revealed that Asp32 or Asp32' in the two forms of mp state similarly exert powerful effects on maintaining both potent and weak inhibitors in the binding pocket of HTLV-1 protease. In the potent inhibitor's unbinding process, His66' was a great supporter that was absent in the weak inhibitor's unbinding pathway. In contrast, in the weak inhibitor's unbinding process, Trp98/Trp98' by pi-pi stacking interactions were unfavorable for the stability of the inhibitor in the binding site. In our opinion, these results will assist in designing more potent and effective inhibitors for the HTLV-1 protease.

All-Atom Simulations of Human ACE2-Spike Protein RBD Complexes for SARS-CoV-2 and Some of its Variants: Nature of Interactions and Free Energy Diagrams for Dissociation of the Protein Complexes.

The spike protein of SARS-CoV-2 is known to interact with the human ACE2 protein via its receptor binding domain (RBD). We have investigated the molecular nature of this interprotein interaction and the associated free energy diagrams for the unbinding of the two proteins for SARS-CoV-2 and some of its known variants through all-atom simulations. The present work involves generation and analysis of 2.5 μs of unbiased and 4.2 μs of biased molecular dynamics trajectories in total for five explicitly solvated RBD-ACE2 systems at full atomic level. First, we have made a comparative analysis of the details of residue-wise specific interactions of the spike protein with ACE2 for SARS-CoV-1 and SARS-CoV-2. It is found that the average numbers of both direct interprotein and water-bridged hydrogen bonds between the RBD and ACE2 are higher for SARS-CoV-2 than SARS-CoV-1. These higher hydrogen bonded interactions are further aided by enhanced nonspecific electrostatic attractions between the two protein surfaces for SARS-CoV-2. The free energy calculations reveal that there is an increase in the free energy barrier by ∼1.5 kcal/mol for the unbinding of RBD from ACE2 for SARS-CoV-2 compared to that for SARS-CoV-1. Subsequently, we considered the RBDs of three variants of SARS-CoV-2, namely N501Y, E484Q/L452R, and N440K. The free energy barrier of protein unbinding for the N501Y variant is found to be ∼4 kcal/mol higher than the wild type SARS-CoV-2 which can be attributed to additional specific interactions involving Tyr501 of RBD and Lys353 and Tyr42 of ACE2 and also enhanced nonspecific electrostatic interaction between the protein surfaces. For the other two mutant variants of E484Q/L452R and N440K, the free energy barrier for protein unbinding increases by ∼2 and ∼1 kcal/mol, respectively, compared with the wild type SARS-CoV-2, which can be attributed to an increase in the number of interprotein hydrogen bonds for the former and also to enhanced positive electrostatic potential on the RBD surfaces for both of the variants. The successive breaking of interprotein hydrogen bonds along the free energy pathway of the unbinding process is also found out for all five systems studied here.

High-sensitive and disposable myocardial infarction biomarker immunosensor with optofluidic microtubule lasing.

The early diagnosis of myocardial infarction can significantly improve the survival rate in emergency treatment, which is mainly implemented by the immunoassay for myocardial infarction biomarkers such as cardiac troponins in blood. In this work, a disposable optofluidic microtubule whispering gallery mode (WGM) immunosensor for label-free cardiac troponin I-C (cTnI-C) complex detection has been proposed and demonstrated with active interrogation enhancement. The disposable microtubule is simply fabricated by a silica capillary with pressurized tapering technology for thin-wall, and the cTnI antibodies are immobilized on the inner wall surface of the microtubule through the self-adherent polydopamine substrate. By configuring the two coupling microfibers, the double-fiber-coupled microtubule cavity can serve as a tunable filter for the mutual-coupled polarimetric fiber ring laser (FRL), whose output laser wavelength is determined by the cTnI-C concentration in the optofluidic microtubule with inherent microfluidic channel. Due to the cyclic-cumulative gain of the FRL, the characteristic resonant peak of optical sensing signal is enhanced in the spectral width compression and the optical signal-to-noise ratio improvement, and therefore the optical immunosensor for cTnI-C can be achieved by tracking the output laser wavelength of the FRL conveniently. The dynamic binding and unbinding process of cTnI-C antigen-antibody is illustrated by monitoring the lasing peak wavelength continuously. Our all-fiber immunosensor demonstrated here has the advantages of fast label-free detection, real-time monitor, high sensitivity and disposable sensing element, which can be an innovative detecting tool in early diagnosis of myocardial infarction.

Jets from main sequence and white dwarf companions during common envelope evolution

ABSTRACT It has long been speculated that jet feedback from accretion on to the companion during a common envelope (CE) event could affect the orbital evolution and envelope unbinding process. We present global 3D hydrodynamical simulations of CE evolution (CEE) that include a jet subgrid model and compare them with an otherwise identical model without a jet. Our binary consists of a 2-M⊙ red giant branch primary and a 1- or 0.5-M⊙ main sequence (MS) or white dwarf (WD) secondary companion modelled as a point particle. We run the simulations for 10 orbits (40 d). Our jet model adds mass at a constant rate $\dot{M}_\mathrm{j}$ of the order of the Eddington rate, with maximum velocity vj of the order of the escape speed, to two spherical sectors with the jet axis perpendicular to the orbital plane. We explore the influence of the jet on orbital evolution, envelope morphology and envelope unbinding, and assess the dependence of the results on the jet mass-loss rate, launch speed, companion mass, opening angle, and accretion rate. In line with our theoretical estimates, jets are choked around the time of first periastron passage and remain choked thereafter. Subsequent to choking, but not before, jets efficiently transfer energy to bound envelope material. This leads to increases in unbound mass of up to $\sim 10{{\ \rm per\ cent}}$, as compared to the simulations without jets. We also estimate the cumulative effects of jets over a full CE phase, finding that jets launched by MS and WD companions are unlikely to dominate envelope unbinding.

Local Ion Densities can Influence Transition Paths of Molecular Binding.

Improper reaction coordinates can pose significant problems for path-based binding free energy calculations. Particularly, omission of long timescale motions can lead to over-estimation of the energetic barriers between the bound and unbound states. Many methods exist to construct the optimal reaction coordinate using a pre-defined basis set of features. Although simulations are typically conducted in explicit solvent, the solvent atoms are often excluded by these feature sets—resulting in little being known about their role in reaction coordinates, and ultimately, their role in determining (un)binding rates and free energies. In this work, analysis is done on an extensive set of host-guest unbinding trajectories, working to characterize differences between high and low probability unbinding trajectories with a focus on solvent-based features, including host-ion interactions, guest-ion interactions and location-dependent ion densities. We find that differences in ion densities as well as guest-ion interactions strongly correlate with differences in the probabilities of reactive paths that are used to determine free energies of (un)binding and play a significant role in the unbinding process.

Assessing models of force-dependent unbinding rates via infrequent metadynamics.

Protein-ligand interactions are crucial for a wide range of physiological processes. Many cellular functions result in these non-covalent "bonds" being mechanically strained, and this can be integral to proper cellular function. Broadly, two classes of force dependence have been observed-slip bonds, where the unbinding rate increases, and catch bonds, where the unbinding rate decreases. Despite much theoretical work, we cannot predict for which protein-ligand pairs, pulling coordinates, and forces a particular rate dependence will appear. Here, we assess the ability of MD simulations combined with enhanced sampling techniques to probe the force dependence of unbinding rates. We show that the infrequent metadynamics technique correctly produces both catch and slip bonding kinetics for model potentials. We then apply it to the well-studied case of a buckyball in a hydrophobic cavity, which appears to exhibit an ideal slip bond. Finally, we compute the force-dependent unbinding rate of biotin-streptavidin. Here, the complex nature of the unbinding process causes the infrequent metadynamics method to begin to break down due to the presence of unbinding intermediates, despite the use of a previously optimized sampling coordinate. Allowing for this limitation, a combination of kinetic and free energy computations predicts an overall slip bond for larger forces consistent with prior experimental results although there are substantial deviations at small forces that require further investigation. This work demonstrates the promise of predicting force-dependent unbinding rates using enhanced sampling MD techniques while also revealing the methodological barriers that must be overcome to tackle more complex targets in the future.

Combined Free-Energy Calculation and Machine Learning Methods for Understanding Ligand Unbinding Kinetics.

The determination of drug residence times, which define the time an inhibitor is in complex with its target, is a fundamental part of the drug discovery process. Synthesis and experimental measurements of kinetic rate constants are, however, expensive and time consuming. In this work, we aimed to obtain drug residence times computationally. Furthermore, we propose a novel algorithm to identify molecular design objectives based on ligand unbinding kinetics. We designed an enhanced sampling technique to accurately predict the free-energy profiles of the ligand unbinding process, focusing on the free-energy barrier for unbinding. Our method first identifies unbinding paths determining a corresponding set of internal coordinates (ICs) that form contacts between the protein and the ligand; it then iteratively updates these interactions during a series of biased molecular dynamics (MD) simulations to reveal the ICs that are important for the whole of the unbinding process. Subsequently, we performed finite-temperature string simulations to obtain the free-energy barrier for unbinding using the set of ICs as a complex reaction coordinate. Importantly, we also aimed to enable the further design of drugs focusing on improved residence times. To this end, we developed a supervised machine learning (ML) approach with inputs from unbiased “downhill” trajectories initiated near the transition state (TS) ensemble of the string unbinding path. We demonstrate that our ML method can identify key ligand–protein interactions driving the system through the TS. Some of the most important drugs for cancer treatment are kinase inhibitors. One of these kinase targets is cyclin-dependent kinase 2 (CDK2), an appealing target for anticancer drug development. Here, we tested our method using two different CDK2 inhibitors for the potential further development of these compounds. We compared the free-energy barriers obtained from our calculations with those observed in available experimental data. We highlighted important interactions at the distal ends of the ligands that can be targeted for improved residence times. Our method provides a new tool to determine unbinding rates and to identify key structural features of the inhibitors that can be used as starting points for novel design strategies in drug discovery.

Ritonavir and xk263 Binding-Unbinding with HIV-1 Protease: Pathways, Energy and Comparison.

Understanding non-covalent biomolecular recognition, which includes drug–protein bound states and their binding/unbinding processes, is of fundamental importance in chemistry, biology, and medicine. Fully revealing the factors that govern the binding/unbinding processes can further assist in designing drugs with desired binding kinetics. HIV protease (HIVp) plays an integral role in the HIV life cycle, so it is a prime target for drug therapy. HIVp has flexible flaps, and the binding pocket can be accessible by a ligand via various pathways. Comparing ligand association and dissociation pathways can help elucidate the ligand–protein interactions such as key residues directly involved in the interaction or specific protein conformations that determine the binding of a ligand under certain pathway(s). Here, we investigated the ligand unbinding process for a slow binder, ritonavir, and a fast binder, xk263, by using unbiased all-atom accelerated molecular dynamics (aMD) simulation with a re-seeding approach and an explicit solvent model. Using ritonavir-HIVp and xk263-HIVp ligand–protein systems as cases, we sampled multiple unbinding pathways for each ligand and observed that the two ligands preferred the same unbinding route. However, ritonavir required a greater HIVp motion to dissociate as compared with xk263, which can leave the binding pocket with little conformational change of HIVp. We also observed that ritonavir unbinding pathways involved residues which are associated with drug resistance and are distal from catalytic site. Analyzing HIVp conformations sampled during both ligand–protein binding and unbinding processes revealed significantly more overlapping HIVp conformations for ritonavir-HIVp rather than xk263-HIVp. However, many HIVp conformations are unique in xk263-HIVp unbinding processes. The findings are consistent with previous findings that xk263 prefers an induced-fit model for binding and unbinding, whereas ritonavir favors a conformation selection model. This study deepens our understanding of the dynamic process of ligand unbinding and provides insights into ligand–protein recognition mechanisms and drug discovery.

Comparative study of the unbinding process of some HTLV-1 protease inhibitors using Unbiased Molecular Dynamics simulation

Comparative study of the unbinding process of some HTLV-1 protease inhibitors using Unbiased Molecular Dynamics simulation