Prediction and analysis of binding affinities for chemically diverse HIV-1 PR inhibitors by the modified SAFE_p approach

One of the biggest challenges in the "in silico" screening of enzyme ligands is to have a protocol that could predict the ligand binding free energies. In our group we have developed a very simple screening function (referred to as solvent accessibility free energy of binding predictor, SAFE_p) which we have applied previously to the study of peptidic HIV-1 protease (HIV-1 PR) inhibitors and later to cyclic urea type HIV-1 PR inhibitors. In this work, we have extended the SAFE_p protocol to a chemically diverse set of HIV-1 PR inhibitors with binding constants that differ by several orders of magnitude. The resulting function is able to reproduce the ranking and in many cases the value of the inhibitor binding affinities for the HIV-1 PR, with accuracy comparable with that of costlier protocols. We also demonstrate that the binding pocket SAFE_p analysis can contribute to the understanding of the physical forces that participate in ligand binding. The analysis tools afforded by our protocol have allowed us to identify an induced fit phenomena mediated by the inhibitor and have demonstrated that larger fragments do not necessarily contribute the most to the binding free energy, an outcome partially brought about by the substantial role the desolvation penalty plays in the energetics of binding. Finally, we have revisited the effect of the Asp dyad protonation state on the predicted binding affinities.     

Does a diol cyclic urea inhibitor of HIV-1 protease bind tighter than its corresponding alcohol form? A study by free energy perturbation and continuum electrostatics calculations

The cyclic urea inhibitors of HIV-1 protease generally have two hydroxyl groups on the seven-membered ring. In this study, free energy perturbation and continuum electrostatic calculations were used to study the contributions of the two hydroxyl groups to the binding affinity and solubility of a cyclic urea inhibitor DMP323. The results indicated that the inhibitor with one hydroxyl group has better binding affinity and solubility than the inhibitor with two hydroxyl groups. Therefore, removal of one hydroxyl group from DMP323 may help to improve the properties of DMP323. This is also likely to be true for other cyclic urea inhibitors. The study also illustrated the difficulty in accurate modeling of the binding affinities of HIV-1 protease inhibitors, which involves many possible protonation states of the two catalytic aspartic acids in the active site of the enzyme.     

An insight to the molecular interactions of the FDA approved HIV PR drugs against L38L↑N↑L PR mutant

The aspartate protease of the human immune deficiency type-1 virus (HIV-1) has become a crucial antiviral target in which many useful antiretroviral inhibitors have been developed. However, it seems the emergence of new HIV-1 PR mutations enhances drug resistance, hence, the available FDA approved drugs show less activity towards the protease. A mutation and insertion designated L38L↑N↑L PR was recently reported from subtype of C-SA HIV-1. An integrated two-layered ONIOM (QM:MM) method was employed in this study to examine the binding affinities of the nine HIV PR inhibitors against this mutant. The computed binding free energies as well as experimental data revealed a reduced inhibitory activity towards the L38L↑N↑L PR in comparison with subtype C-SA HIV-1 PR. This observation suggests that the insertion and mutations significantly affect the binding affinities or characteristics of the HIV PIs and/or parent PR. The same trend for the computational binding free energies was observed for eight of the nine inhibitors with respect to the experimental binding free energies. The outcome of this study shows that ONIOM method can be used as a reliable computational approach to rationalize lead compounds against specific targets. The nature of the intermolecular interactions in terms of the host–guest hydrogen bond interactions is discussed using the atoms in molecules (AIM) analysis. Natural bond orbital analysis was also used to determine the extent of charge transfer between the QM region of the L38L↑N↑L PR enzyme and FDA approved drugs. AIM analysis showed that the interaction between the QM region of the L38L↑N↑L PR and FDA approved drugs are electrostatic dominant, the bond stability computed from the NBO analysis supports the results from the AIM application. Future studies will focus on the improvement of the computational model by considering explicit water molecules in the active pocket. We believe that this approach has the potential to provide information that will aid in the design of much improved HIV-1 PR antiviral drugs.     

Automated molecular simulation based binding affinity calculator for ligand-bound HIV-1 proteases

The successful application of high throughput molecular simulations to determine biochemical properties would be of great importance to the biomedical community if such simulations could be turned around in a clinically relevant timescale. An important example is the determination of antiretroviral inhibitor efficacy against varying strains of HIV through calculation of drug-protein binding affinities. We describe the Binding Affinity Calculator (BAC), a tool for the automated calculation of HIV-1 protease-ligand binding affinities. The tool employs fully atomistic molecular simulations alongside the well established molecular mechanics Poisson-Boltzmann solvent accessible surface area (MMPBSA) free energy methodology to enable the calculation of the binding free energy of several ligand-protease complexes, including all nine FDA approved inhibitors of HIV-1 protease and seven of the natural substrates cleaved by the protease. This enables the efficacy of these inhibitors to be ranked across several mutant strains of the protease relative to the wildtype. BAC is a tool that utilizes the power provided by a computational grid to automate all of the stages required to compute free energies of binding: model preparation, equilibration, simulation, postprocessing, and data-marshaling around the generally widely distributed compute resources utilized. Such automation enables the molecular dynamics methodology to be used in a high throughput manner not achievable by manual methods. This paper describes the architecture and workflow management of BAC and the function of each of its components. Given adequate compute resources, BAC can yield quantitative information regarding drug resistance at the molecular level within 96 h. Such a timescale is of direct clinical relevance and can assist in decision support for the assessment of patient-specific optimal drug treatment and the subsequent response to therapy for any given genotype.     

Energetic basis for drug resistance of HIV-1 protease mutants against amprenavir

Amprenavir (APV) is a high affinity (0.15 nM) HIV-1 protease (PR) inhibitor. However, the affinities of the drug resistant protease variants V32I, I50V, I54V, I54M, I84V and L90M to amprenavir are decreased 3 to 30-fold compared to the wild-type. In this work, the popular molecular mechanics Poisson-Boltzmann surface area method has been used to investigate the effectiveness of amprenavir against the wild-type and these mutated protease variants. Our results reveal that the protonation state of Asp25/Asp25′ strongly affects the dynamics, the overall affinity and the interactions of the inhibitor with individual residues. We emphasize that, in contrast to what is often assumed, the protonation state may not be inferred from the affinities but requires pK_a calculations. At neutral pH, Asp25 and Asp25′ are ionized or protonated, respectively, as suggested from pK_a calculations. This protonation state was thus mainly considered in our study. Mutation induced changes in binding affinities are in agreement with the experimental findings. The decomposition of the binding free energy reveals the mechanisms underlying binding and drug resistance. Drug resistance arises from an increase in the energetic contribution from the van der Waals interactions between APV and PR (V32I, I50V, and I84V mutant) or a rise in the energetic contribution from the electrostatic interactions between the inhibitor and its target (I54M and I54V mutant). For the V32I mutant, also an increased free energy for the polar solvation contributes to the drug resistance. For the L90M mutant, a rise in the van der Waals energy for APV-PR interactions is compensated by a decrease in the polar solvation free energy such that the net binding affinity remains unchanged. Detailed understanding of the molecular forces governing binding and drug resistance might assist in the design of new inhibitors against HIV-1 PR variants that are resistant against current drugs.     

新型二氟甲基磷酸类酪氨酸蛋白磷酸酯酶1B抑制剂的分子动力学模拟和结合自由能计算

通过分子对接建立了一系列含二氟甲基磷酸基团(DFMP)或二氟甲基硫酸基团(DFMS)的抑制剂与酪氨酸蛋白磷酸酯酶1B(PTP1B)的相互作用模式, 并通过1 ns的分子动力学模拟和molecular mechanics/generalized Born surface area (MM/GBSA)方法计算了其结合自由能. 计算获得的结合自由能排序和抑制剂与靶酶间结合能力排序一致; 通过基于主方程的自由能计算方法, 获得了抑制剂与靶酶残基间相互作用的信息, 这些信息显示DFMP/DFMS基团的负电荷中心与PTP1B的221位精氨酸正电荷中心之间的静电相互作用强弱决定了此类抑制剂的活性, 进一步的分析还显示位于DFMP/DFMS基团中的氟原子或其他具有适当原子半径的氢键供体原子会增进此类抑制剂与PTP1B活性位点的结合能力.     

The Normal-Mode Entropy in the MM/GBSA Method: Effect of System Truncation, Buffer Region, and Dielectric Constant

We have performed a systematic study of the entropy term in the MM/GBSA (molecular mechanics combined with generalized Born and surface-area solvation) approach to calculate ligand-binding affinities. The entropies are calculated by a normal-mode analysis of harmonic frequencies from minimized snapshots of molecular dynamics simulations. For computational reasons, these calculations have normally been performed on truncated systems. We have studied the binding of eight inhibitors of blood clotting factor Xa, nine ligands of ferritin, and two ligands of HIV-1 protease and show that removing protein residues with distances larger than 8-16 ? to the ligand, including a 4 ? shell of fixed protein residues and water molecules, change the absolute entropies by 1-5 kJ/mol on average. However, the change is systematic, so relative entropies for different ligands change by only 0.7-1.6 kJ/mol on average. Consequently, entropies from truncated systems give relative binding affinities that are identical to those obtained for the whole protein within statistical uncertainty (1-2 kJ/mol). We have also tested to use a distance-dependent dielectric constant in the minimization and frequency calculation (ε = 4r), but it typically gives slightly different entropies and poorer binding affinities. Therefore, we recommend entropies calculated with the smallest truncation radius (8 ?) and ε =1. Such an approach also gives an improved precision for the calculated binding free energies.     

分子动力学模拟研究质子化态在HIV-1 Protease-Indinavir复合物中的作用

I 型人体免疫缺陷病毒(HIV-1)蛋白酶中Asp25/Asp25'的质子化对于理论研究HIV-1 蛋白酶和抑制剂的作用机制
以及氨基酸变异对抗药性的影响有重要意义. 分别对Protease-Indinavir (PR-IDV)复合物的六种可能的质子化态进行了
5 ns 的分子动力学模拟, 分析了不同状态对动力学特征和结构的影响, 用molecular mechanics/Possion-Boltzman surface
area (MM-PBSA)方法计算了PR 和IDV 在各种状态下的结合自由能. 计算结果说明A 链Asp25 的OD2 的质子化是最
为可能的状态. 对PR-IDV 复合物中起到媒介作用的水分子与PR-IDV 复合物形成的氢键进行了分析, 分析结果说明不
同的质子化态对水分子在PR-IDV 复合物中所起的媒介作用没有影响, 这一结果与我们先前对PR-BEA369 复合物的研
究不同. 我们的研究结果为更高效的PR 抑制剂的设计以及PR 氨基酸变异对药物抗药性的研究提供了理论上的指导.     

Continuum solvation models in the linear interaction energy method

The linear interaction energy (LIE) method in combination with two different continuum solvent models has been applied to calculate protein-ligand binding free energies for a set of inhibitors against the malarial aspartic protease plasmepsin II. Ligand-water interaction energies are calculated from both Poisson-Boltzmann (PB) and Generalized Born (GB) continuum models using snapshots from explicit solvent simulations of the ligand and protein-ligand complex. These are compared to explicit solvent calculations, and we find close agreement between the explicit water and PB solvation models. The GB model overestimates the change in solvation energy, and this is caused by consistent underestimation of the effective Born radii in the protein-ligand complex. The explicit solvent LIE calculations and LIE-PB, with our standard parametrization, reproduce absolute experimental binding free energies with an average unsigned error of 0.5 and 0.7 kcal/mol, respectively. The LIE-GB method, however, requires a constant offset to approach the same level of accuracy.     

分子动力学研究V82A和L90M变异对HIV PR-IDV复合物的影响

为了说明V82A和L90M变异对蛋白酶(PR)和茚地那韦(IDV)复合物的影响,进行了5.5ns的MD模拟.用MM-PBSA方法计算了体系的结合自由能,计算和实验结果一致.分解自由能为不同能量项说明,这两个变异引起熵的贡献变化大于焓的贡献变化.分解自由能到每个残基说明Wild,V82A和L90M具有相似的结合模式,结合能的贡献主要来源于A28/A28',I50/I50'和I84/I84'这六个残基组,详细分析了Wild和IDV的结合模式,对比分析了V82A和L90M变异引起结合模式的细小变化.V82A变异引起结合模式的变化是由于变异后位阻减小导致的.L90M变异引起D25和L90间的作用增强并引起结合模式的细小变化.研究结果有助于更好地理解变异对抑制剂和HIV-1PR结合模式的影响,并可以用来帮助设计更高效的PR抑制剂.     

Characterizing the protonation states of the catalytic residues in apo and substrate-bound human T-cell leukemia virus type 1 protease

Human T-cell leukemia virus type 1 (HTLV-1) protease is an attractive target when developing inhibitors to treat HTLV-1 associated diseases. To study the catalytic mechanism and design novel HTLV-1 protease inhibitors, the protonation states of the two catalytic aspartic acid residues must be determined. Free energy simulations have been conducted to study the proton transfer reaction between the catalytic residues of HTLV-1 protease using a combined quantum mechanical and molecular mechanical (QM/MM) molecular dynamics simulation. The free energy profiles for the reaction in the apo-enzyme and in an enzyme – substrate complex have been obtained. In the apo-enzyme, the two catalytic residues are chemically equivalent and are expected to be both unprotonated. Upon substrate binding, the catalytic residues of HTLV-1 protease evolve to a singly protonated state, in which the OD1 of Asp32 is protonated and forms a hydrogen bond with the OD1 of Asp32′, which is unprotonated. The HTLV-1 protease–substrate complex structure obtained from this simulation can serve as the Michaelis complex structure for further mechanistic studies of HTLV-1 protease while providing a receptor structure with the correct protonation states for the active site residues toward the design of novel HTLV-1 protease inhibitors through virtual screening.     

胰蛋白酶和苯酰氨类抑制剂结合自由能的预测

用基于线性响应近似的自由能预测方法计算胰蛋白酶和苯酰氨类抑制剂的结合 自由能。计算结果表明,单参数,双参数和三参数模型具有相似的线性回归系数, 但三参数和双参数模型的交互验证回归系数要明显优于单参数模型。从预测能力来 看,双参数模型和三参数模型都能够很好地预测测试集中抑制剂的结合自由能,其 中双参数模型预测的结果要略优于三参数模型的预测结果。对测试集中的抑制剂, 双参数模型预测得到的预测自由能和实际自由能之间平均绝对误差仅为1.15 kJ/mol。自由能计算模型以及分子动力学轨迹能很好地解释抑制剂结构和活性的 关系,为药物设计提供了重要的结构信息。     

A combined treatment of hydration and dynamical effects for the modeling of host–guest binding thermodynamics: the SAMPL5 blinded challenge

As part of the SAMPL5 blinded experiment, we computed the absolute binding free energies of 22 host–guest complexes employing a novel approach based on the BEDAM single-decoupling alchemical free energy protocol with parallel replica exchange conformational sampling and the AGBNP2 implicit solvation model specifically customized to treat the effect of water displacement as modeled by the Hydration Site Analysis method with explicit solvation. Initial predictions were affected by the lack of treatment of ionic charge screening, which is very significant for these highly charged hosts, and resulted in poor relative ranking of negatively versus positively charged guests. Binding free energies obtained with Debye–Hückel treatment of salt effects were in good agreement with experimental measurements. Water displacement effects contributed favorably and very significantly to the observed binding affinities; without it, the modeling predictions would have grossly underestimated binding. The work validates the implicit/explicit solvation approach employed here and it shows that comprehensive physical models can be effective at predicting binding affinities of molecular complexes requiring accurate treatment of conformational dynamics and hydration.     

Effect of atomic charge,solvation, entropy,and ligand protonation state on MM-PB(GB)SA binding energies of HIV protease

The molecular mechanics-Poisson-Boltzmann surface area (MM-PBSA) and MM-generalized-Born surface area (MM-GBSA) approaches are commonly used in molecular modeling and drug design. Four critical aspects of these approaches have been investigated for their effect on calculated binding energies: (1) the atomic partial charge method used to parameterize the ligand force field, (2) the method used to calculate the solvation free energy, (3) inclusion of entropy estimates, and (4) the protonation state of the ligand. HIV protease has been used as a test case with six structurally different inhibitors covering a broad range of binding strength to assess the effect of these four parameters. Atomic charge methods are demonstrated to effect both the molecular dynamics (MD) simulation and MM-PB(GB)SA binding energy calculation, with a greater effect on the MD simulation. Coefficients of determination and Spearman rank coefficients were used to quantify the performance of the MM-PB(GB)SA methods relative to the experimental data. In general, better performance was achieved using (i) atomic charge models that produced smaller mean absolute atomic charges (Gasteiger, HF/STO-3G and B3LYP/cc-pVTZ), (ii) the MM-GBSA approach over MM-PBSA, while (iii) inclusion of entropy had a slightly positive effect on correlations with experiment. Accurate representation of the ligand protonation state was found to be important. It is demonstrated that these approaches can distinguish ligands according to binding strength, underlining the usefulness of these approaches in computer-aided drug design. © 2012 Wiley Periodicals, Inc.     

Binding of novel fullerene inhibitors to HIV-1 protease: insight through molecular dynamics and molecular mechanics Poisson–Boltzmann surface area calculations

The objectives of this study include the design of a series of novel fullerene-based inhibitors for HIV-1 protease (HIV-1 PR), by employing two strategies that can also be applied to the design of inhibitors for any other target. Additionally, the interactions which contribute to the observed exceptionally high binding free energies were analyzed. In particular, we investigated: (1) hydrogen bonding (H-bond) interactions between specific fullerene derivatives and the protease, (2) the regions of HIV-1 PR that play a significant role in binding, (3) protease changes upon binding and (4) various contributions to the binding free energy, in order to identify the most significant of them. This study has been performed by employing a docking technique, two 3D-QSAR models, molecular dynamics (MD) simulations and the molecular mechanics Poisson–Boltzmann surface area (MM–PBSA) method. Our computed binding free energies are in satisfactory agreement with the experimental results. The suitability of specific fullerene derivatives as drug candidates was further enhanced, after ADMET (absorption, distribution, metabolism, excretion and toxicity) properties have been estimated to be promising. The outcomes of this study revealed important protein–ligand interaction patterns that may lead towards the development of novel, potent HIV-1 PR inhibitors.     

Rapid and accurate prediction of binding free energies for saquinavir-bound HIV-1 proteases

To explain drug resistance by computer simulations at the molecular level, we first have to assess the accuracy of theoretical predictions. Herein we report an application of the molecular mechanics Poisson-Boltzmann surface area (MM/PBSA) technique to the ranking of binding affinities of the inhibitor saquinavir with the wild type (WT) and three resistant mutants of HIV-1 protease: L90M, G48V, and G48V/L90M. For each ligand-protein complex we report 10 ns of fully unrestrained molecular dynamics (MD) simulations with explicit solvent. We investigate convergence, internal consistency, and model dependency of MM/PBSA ligand binding energies. Converged enthalpy and entropy estimates produce ligand binding affinities within 1.5 kcal/mol of experimental values, with a remarkable level of correlation to the experimentally observed ranking of resistance levels. A detailed analysis of the enthalpic/entropic balance of drug-protease interactions explains resistance in L90M in terms of a higher vibrational entropy than in the WT complex, while G48V disrupts critical hydrogen bonds at the inhibitor's binding site and produces an altered, more unfavorable balance of Coulomb and polar desolvation energies.     

Molecular dynamics simulation and free energy calculation studies of the binding mechanism of allosteric inhibitors with p38α MAP kinase

p38 MAP kinase is a promising target for anti-inflammatory treatment. The classical kinase inhibitors imatinib and sorafenib as well as BI-1 and BIRB-796 were reported to bind in the DFG-out form of human p38α, known as type II or allosteric kinase inhibitors. Although DFG-out conformation has attracted great interest in the design of type II kinase inhibitors, the structural requirements for binding and mechanism of stabilization of DFG-out conformation remain unclear. As allosteric inhibition is important to the selectivity of kinase inhibitor, herein the binding modes of imatinib, sorafenib, BI-1 and BIRB-796 to p38α were investigated by molecular dynamics simulation. Binding free energies were calculated by molecular mechanics/Poisson-Boltzmann surface area method. The predicted binding affinities can give a good explanation of the activity difference of the studied inhibitors. Furthermore, binding free energies decomposition analysis and further structural analysis indicate that the dominating effect of van der Waals interaction drives the binding process, and key residues, such as Lys53, Gly71, Leu75, Ile84, Thr106, Met109, Leu167, Asp168, and Phe169, play important roles by forming hydrogen bond, salt bridge, and hydrophobic interactions with the DFG-out conformation of p38α. Finally, we also conducted a detailed analysis of BI-1, imatinib, and sorafenib binding to p38α in comparison with BIRB-796 exploited for gaining potency as well as selectivity of p38 inhibitors. These results are expected to be useful for future rational design of novel type II p38 inhibitors.     

Calculation of solvation and binding free energy differences between VX-478 and its analogs by free energy perturbation and AMSOL methods

Summary VX-478 belongs to a novel class of HIV-1 protease inhibitors that are based on N,N-disubstituted benzene sulfonamides. Force field parameters for the N,N-dialkyl benzene sulfonamide moiety have been assembled from the literature and from our own ab initio calculations. These parameters were employed to calculate solvation and binding free energy differences between VX-478 and two analogs. The free energy perturbation method has been used to determine these differences using two approaches. In the first approach, intergroup interaction terms only were included in the calculation of free energies (as in most reports of free energy calculations using AMBER). In the second approach, both the inter- and intragroup interaction terms were included. The results obtained with the two approaches are in excellent agreement with each other and are also in close agreement with the experimental results. The solvation free energies of N,N-dimethyl benzene sulfonamide derivatives (truncated models of the inhibitors), calculated using continuum solvation (AMSOL) methods, are found to be in qualitative agreement with the experimental and free energy perturbation results. The binding and solvation free energy results are discussed in the context of structure-based drug design to show how physicochemical properties (for example aqueous solubilities and bioavailabilities) of these HIV-1 protease inhibitors were improved, while maintaining their inhibitory potency.     

Binding free energy contributions of interfacial waters in HIV-1 protease/inhibitor complexes

Water molecules are commonly observed in crystal structures of protein-ligand complexes where they mediate protein-ligand binding. It is of considerable theoretical and practical importance to determine quantitatively the individual free energy contributions of these interfacial water molecules to protein-ligand binding and to elucidate factors that influence them. The double-decoupling free energy molecular dynamics simulation method has been used to calculate the binding free energy contribution for each of the four interfacial water molecules observed in the crystal structure of HIV-1 protease complexed with KNI-272, a potent inhibitor. While two of these water molecules contribute significantly to the binding free energy, the other two have close to zero contribution. It was further observed that the protonation states of two catalytic aspartate residues, Asp25 and Asp125, strongly influence the free energy contribution of a conserved water molecule Wat301 and that different inhibitors significantly influence the free energy contribution of Wat301. Our results have important implications on our understanding of the role of interfacial water molecules in protein-ligand binding and to structure-based drug design aimed at incorporating these interfacial water molecules into ligands.