How Good is Jarzynski's Equality for Computer-Aided Drug Design?

Multisite λ-Dynamics for Protein-DNA Binding Affinity Prediction

Transcription factors (TFs) regulate gene expression by binding to specific DNA sequences, playing critical roles in cellular processes and disease pathways. Computational methods, particularly λ-Dynamics, offer a promising approach for predicting TF relative binding affinities. This study evaluates the effectiveness of different λ-Dynamics perturbation schemes in determining binding free energy changes (ΔΔGb) of the WRKY transcription factor upon mutating its W-box binding site (GGTCAA) to a nonspecific sequence (GATAAA). Among the schemes tested, the single λ per base pair protocol demonstrated the fastest convergence and highest precision. Extending this protocol to additional mutants (GGTCCG and GGACAA) yielded ΔΔGb values that successfully ranked binding affinities, showcasing its strong potential for high-throughput screening of DNA binding sites.

Integrated computational biophysics approach for drug discovery against Nipah virus

The Nipah virus (NiV) poses a pressing global threat to public health due to its high mortality rate, multiple modes of transmission, and lack of effective treatments. NiV glycoprotein G (NiV-G) emerges as a promising target for the discovery of NiV drugs because of its essential role in viral entry and membrane fusion. Therefore, in this study, we applied an integrated computational and biophysics approach to identify potential inhibitors of NiV-G within a curated dataset of Peruvian phytochemicals. The virtual screening results indicated that these compounds could represent a natural source of potential NiV-G inhibitors with ΔG values ranging from -8 to -11 kcal/mol. Among them, procyanidin B2, B3, B7, and C1 exhibited the highest binding affinities and formed the most molecular interactions with NiV-G. Molecular dynamics simulations revealed the induced-fit mechanism of NiV-G pocket interaction with these procyanidins, primarily driven by its hydrophobic nature. Non-equilibrium free energy calculations were used to determine binding affinities, highlighting Procyanidin B3 and B2 as the ligands with the most substantial interactions. In general, this work underscores the potential of Peruvian phytochemicals, particularly procyanidins B2, B3, B7, and C1, as lead compounds for developing anti-NiV drugs through an integrated computational biophysics approach.

Treatment of flexibility of protein backbone in simulations of protein-ligand interactions using steered molecular dynamics

To ensure that an external force can break the interaction between a protein and a ligand, the steered molecular dynamics simulation requires a harmonic restrained potential applied to the protein backbone. A usual practice is that all or a certain number of protein's heavy atoms or Cα atoms are fixed, being restrained by a small force. This present study reveals that while fixing both either all heavy atoms and or all Cα atoms is not a good approach, while fixing a too small number of few atoms sometimes cannot prevent the protein from rotating under the influence of the bulk water layer, and the pulled molecule may smack into the wall of the active site. We found that restraining the Cα atoms under certain conditions is more relevant. Thus, we would propose an alternative solution in which only the Cα atoms of the protein at a distance larger than 1.2 nm from the ligand are restrained. A more flexible, but not too flexible, protein will be expected to lead to a more natural release of the ligand.

Crystallographic and Computational Insights into Isoform-Selective Dynamics in Nitric Oxide Synthase

In our efforts to develop inhibitors selective for neuronal nitric oxide synthase (nNOS) over endothelial nitric oxide synthase (eNOS), we found that nNOS can undergo conformational changes in response to inhibitor binding that does not readily occur in eNOS. One change involves movement of a conserved tyrosine, which hydrogen bonds to one of the heme propionates, but in the presence of an inhibitor, changes conformation, enabling part of the inhibitor to hydrogen bond with the heme propionate. This movement does not occur as readily in eNOS and may account for the reason why these inhibitors bind more tightly to nNOS. A second structural change occurs upon the binding of a second inhibitor molecule to nNOS, displacing the pterin cofactor. Binding of this second site inhibitor requires structural changes at the dimer interface, which also occurs more readily in nNOS than in eNOS. Here, we used a combination of crystallography, mutagenesis, and computational methods to better understand the structural basis for these differences in NOS inhibitor binding. Computational results show that a conserved tyrosine near the primary inhibitor binding site is anchored more tightly in eNOS than in nNOS, allowing for less flexibility of this residue. We also find that the inefficiency of eNOS to bind a second inhibitor molecule is likely due to the tighter dimer interface in eNOS compared with nNOS. This study provides a better understanding of how subtle structural differences in NOS isoforms can result in substantial dynamic differences that can be exploited in the development of isoform-selective inhibitors.

Fluctuation Relations to Calculate Protein Redox Potentials from Molecular Dynamics Simulations

The tunable design of protein redox potentials promises to open a range of applications in biotechnology and catalysis. Here, we introduce a method to calculate redox potential changes by combining fluctuation relations with molecular dynamics simulations. It involves the simulation of reduced and oxidized states, followed by the instantaneous conversion between them. Energy differences introduced by the perturbations are obtained using the Kubo-Onsager approach. Using a detailed fluctuation relation coupled with Bayesian inference, these are postprocessed into estimates for the redox potentials in an efficient manner. This new method, denoted MD + CB, is tested on a de novo four-helix bundle heme protein (the m4D2 "maquette") and five designed mutants, including some mutants characterized experimentally in this work. The MD + CB approach is found to perform reliably, giving redox potential shifts with reasonably good correlation (0.85) to the experimental values for the mutants. The MD + CB approach also compares well with redox potential shift predictions using a continuum electrostatic method. The estimation method employed within the MD + CB approach is straightforwardly transferable to standard equilibrium MD simulations and holds promise for redox protein engineering and design applications.

Nonspecific vs. specific DNA binding free energetics of a transcription factor domain protein

Transcription factor (TF) proteins regulate gene expression by binding to specific sites on the genome. In the facilitated diffusion model, an optimized search process is achieved by the TF alternating between 3D diffusion in the bulk and 1D diffusion along DNA. While undergoing 1D diffusion, the protein can switch from a search mode for fast diffusion along nonspecific DNA to a recognition mode for stable binding to specific DNA. It was recently noticed that, for a small TF domain protein, reorientations on DNA happen between the nonspecific and specific DNA binding. We here conducted all-atom molecular dynamics simulations with steering forces to reveal the protein-DNA binding free energetics, confirming that the search and recognition modes are distinguished primarily by protein orientations on the DNA. As the binding free energy difference between the specific and nonspecific DNA system slightly deviates from that being estimated directly from dissociation constants on 15-bp DNA constructs, we hypothesize that the discrepancy can come from DNA sequences flanking the 6-bp central binding sites that impact on the dissociation kinetics measurements. The hypothesis is supported by a simplified spherical protein-DNA model along with stochastic simulations and kinetic modeling.

Dissociation Rate Calculation via Constant-Force Steered Molecular Dynamics Simulation

Steered molecular dynamics (SMD) simulations are used to study molecular dissociation events by applying a harmonic force to the molecules and pulling them at a constant velocity. Instead of constant-velocity pulling, we use a constant force: the constant-force SMD (CF-SMD) simulation. The CF-SMD simulation employs a constant force to reduce the activation barrier of molecular dissociation, thereby enhancing the dissociation event. Here, we present the capability of the CF-SMD simulation to estimate the dissociation time at equilibrium. We performed all-atom CF-SMD simulations for NaCl and protein-ligand systems, producing dissociation time at various forces. We extrapolated these values to the dissociation rate without a constant force using Bell's model or the Dudko-Hummer-Szabo model. We demonstrate that the CF-SMD simulations with the models predicted the dissociation time in equilibrium. A CF-SMD simulation is a powerful tool for estimating the dissociation rate in a direct and computationally efficient manner.

Improving the Potential of Mean Force and Nonequilibrium Pulling Simulations by Simultaneous Alchemical Modifications

We present an approach combining alchemical modifications and physical-pathway methods to calculate absolute binding free energies. The employed physical-pathway method is either a stratified umbrella sampling to calculate a potential of mean force or nonequilibrium pulling. We devised two basic approaches: the simultaneous approach (S-approach), where, along the physical unbinding pathway, an alchemical transformation of ligand-protein interactions is installed and deinstalled, and the prior-plus-simultaneous approach (PPS-approach), where, prior to the physical-pathway simulation, an alchemical transformation of ligand-protein interactions is installed in the binding site and deinstalled during the physical-pathway simulation. Using a mutant of T4 lysozyme with a benzene ligand as an example, we show that installation and deinstallation of soft-core interactions concurrent with physical ligand unbinding (S-approach) allow successful potential of mean force calculations and nonequilibrium pulling simulations despite the problems posed by the occluded nature of the lysozyme binding pocket. Good agreement between the potential of the mean-force-based S-approach and double decoupling simulations as well as a remarkable efficiency and accuracy of the nonequilibrium-pulling-based S-approach is found. The latter turned out to be more compute-efficient than the potential of mean force calculation by approximately 70%. Furthermore, we illustrate the merits of reducing ligand-protein interactions prior to potential of mean force calculations using the murine double minute homologue protein MDM2 with a p53-derived peptide ligand (PPS-approach). Here, the problem of breaking strong interactions in the binding pocket is transferred to a prior alchemical transformation that reduces the free-energy barrier between the bound and unbound state in the potential of mean force. Besides, disentangling physical ligand displacement from the deinstallation of ligand-protein interactions was seen to allow a more uniform sampling of distance histograms in the umbrella sampling. In the future, physical ligand unbinding combined with simultaneous alchemical modifications may prove useful in the calculation of protein-protein binding free energies, where sampling problems posed by multiple, possibly sticky interactions and potential steric clashes can thus be reduced.

Structural analysis of P450 AmphL from Streptomyces nodosus provides insights into substrate selectivity of polyene macrolide antibiotic biosynthetic P450s

AmphL is a cytochrome P450 enzyme that catalyzes the C8 oxidation of 8-deoxyamphotericin B to the polyene macrolide antibiotic, amphotericin B. To understand this substrate selectivity, we solved the crystal structure of AmphL to a resolution of 2.0 Å in complex with amphotericin B and performed molecular dynamics (MD) simulations. A detailed comparison with the closely related P450, PimD, which catalyzes the epoxidation of 4,5-desepoxypimaricin to the macrolide antibiotic, pimaricin, reveals key catalytic structural features responsible for stereo- and regio-selective oxidation. Both P450s have a similar access channel that runs parallel to the active site I helix over the surface of the heme. Molecular dynamics simulations of substrate binding reveal PimD can “pull” substrates further into the P450 access channel owing to additional electrostatic interactions between the protein and the carboxyl group attached to the hemiketal ring of 4,5-desepoxypimaricin. This substrate interaction is absent in AmphL although the additional substrate -OH groups in 8-deoxyamphotericin B help to correctly position the substrate for C8 oxidation. Simulations of the oxy-complex indicates that these -OH groups may also participate in a proton relay network required for O₂ activation as has been suggested for two other macrolide P450s, PimD and P450eryF. These findings provide experimentally testable models that can potentially contribute to a new generation of novel macrolide antibiotics with enhanced antifungal and/or antiprotozoal efficacy.

Dependence of Work on the Pulling Speed in Mechanical Ligand Unbinding

In single-molecule force spectroscopy, the rupture force F_max required for mechanical unfolding of a biomolecule or for pulling a ligand out of a binding site depends on the pulling speed V and, in the linear Bell-Evans regime, F_max ∼ ln(V). Recently, it has been found that non-equilibrium work W is better than F_max in describing relative ligand binding affinity, but the dependence of W on V remains unknown. In this paper, we developed an analytical theory showing that in the linear regime, W ∼ c₁ ln(V) + c₂ ln²(V), where c₁ and c₂ are constants. This quadratic dependence was also confirmed by all-atom steered molecular dynamics simulations of protein-ligand complexes. Although our theory was developed for ligand unbinding, it is also applicable to other processes, such as mechanical unfolding of proteins and other biomolecules, due to its universality.

Properties of DNA- and Protein-Scaffolded Lipid Nanodiscs

The properties of natural lipid bilayers are vital to the regulation of many membrane proteins. Scaffolded nanodiscs provide an in vitro lipid bilayer platform to host membrane proteins in an environment that approximates native lipid bilayers. However, the properties of scaffold-enclosed bilayers may depart significantly from those of bulk cellular membranes. Therefore, to improve the usefulness of nanodiscs it is essential to understand the properties of lipids restricted by scaffolds. We used computational molecular dynamics and modeling approaches to understand the effects of nanodisc size, scaffold type (DNA or protein), and hydrophobic modification of DNA scaffolds on bilayer stability and degree to which the properties of enclosed bilayers approximate bulk bilayers. With respect to achieving bulk bilayer behavior, we found that charge neutralization of DNA scaffolds was more important than the total hydrophobic content of their modifications: bilayer properties were better for scaffolds having a large number of short alkyl chains than those having fewer long alkyl chains. Further, complete charge neutralization of DNA scaffolds enabled better lipid binding, and more stable bilayers, as shown by steered molecular dynamics simulations that measured the force required to dislodge scaffolds from lipid bilayer patches. Considered together, our simulations provide a guide to the design of DNA-scaffolded nanodiscs suitable for studying membrane proteins.

Structure-based lead optimization of herbal medicine rutin for inhibiting SARS-CoV-2's main protease

Coronavirus disease 2019 (COVID-19) is an ongoing global pandemic with very limited specific treatments. To fight COVID-19, various traditional antiviral medicines have been prescribed in China to infected patients with mild to moderate symptoms and received unexpected success in controlling the disease. However, the molecular mechanisms of how these herbal medicines interact with the SARS-CoV-2 virus that causes COVID-19 have remained elusive. It is well known that the main protease (Mpro) of SARS-CoV-2 plays an important role in maturation of many viral proteins such as the RNA-dependent RNA polymerase. Here, we explore the underlying molecular mechanisms of the computationally determined top candidate, namely, rutin which is a key component in many traditional antiviral medicines such as Lianhuaqinwen and Shuanghuanlian, for inhibiting the viral target-Mpro. Using in silico methods (docking and molecular dynamics simulations), we revealed the dynamics and energetics of rutin when interacting with the Mpro of SARS-CoV-2, suggesting that the highly hydrophilic rutin molecule can be bound inside the Mpro's pocket (active site) and possibly inhibit its biological functions. In addition, we optimized the structure of rutin and designed two more hydrophobic analogs, M1 and M2, which satisfy the rule of five for western medicines and demonstrated that they (M2 in particular) possess much stronger binding affinities to the SARS-COV-2s Mpro than rutin, due to the enhanced hydrophobic interaction as well as more hydrogen bonds. Therefore, our results provide invaluable insights into the mechanism of a ligand's binding inside the Mpro and shed light on future structure-based designs of high-potent inhibitors for SARS-CoV-2 Mpro.