分子动力学模拟研究质子化态在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 氨基酸变异对药物抗药性的研究提供了理论上的指导.     

Z_(Aβ3)和Aβ_(16–40)亲和作用的分子机理解析

淀粉样多肽(amyloid-βpeptide,Aβ)聚集是引起阿尔兹海默症(Alzheimer's disease,AD)的主要原因。开发Aβ聚集抑制剂是治疗AD的最有效手段之一。利用噬菌体展示技术筛选出来的Z_(Aβ3)蛋白质能够有效抑制Aβ聚集,但Z_(Aβ3)和Aβ之间的作用区域和关键氨基酸残基尚不清楚。针对此问题,本研究利用分子动力学模拟、MM-PBSA自由能计算和分解方法研究了Z_(Aβ3)-Aβ_(16–40)复合物之间的相互作用机制。结果表明,Z_(Aβ3)的β-股和Aβ_(16–40)之间的亲和作用占主导,而Z_(Aβ3)的α-螺旋贡献很小。利用分子力学-帕松波尔茨曼溶剂可及化表面积方法(MM-PBSA)自由能分解发现Z_(Aβ3)的热点残基为E15、I16、V17、Y18、L19、P20、N21和L22,而Aβ_(16–40)的热点残基为F19、F20、A21、E22、D23、K28、I31、I32、G33、L34、M35、V36、G38和V40。Z_(Aβ3)通过将发夹型Aβ单体包埋在α-螺旋围成的疏水性腔体内来阻碍Aβ聚集。这种结合模式为设计高效的Aβ蛋白质类抑制剂提供了三个基本要素:高亲和性的结合片段(β-股)、附属结构(α-螺旋)和通过二硫键形成的稳定构象。高亲和性结合片段能竞争性地与Aβ单体结合,附属结构α-螺旋可以阻碍其它Aβ单体靠近,而稳定的构象是上述两种要素发挥作用的基础,三者协同作用可以有效地抑制Aβ聚集。     

乙酰胆碱酯酶AChE与1,7-二氮杂咔唑抑制剂的作用机理的分子动力学模拟

通过分子对接、分子动力学(MD)模拟以及成键自由能分析方法,从原子水平上模拟研究了3种1,7-二氮杂咔唑衍生物(分别记为M1、M2和M3)与ACh E的结合模式及相互作用机理,分析和讨论了研究体系的静电相互作用和范德华相互作用(vd W)。用MM-PBSA方法计算的3种抑制剂与ACh E之间的结合自由能与抑制剂的实验生物活性数据(IC50值)相对应。分析结果表明,残基S286与抑制剂之间形成的氢键作用有利于抑制剂与ACh E之间的结合。范德华相互作用,尤其是抑制剂与关键残基W279和Y334的作用,对抑制剂与ACh E之间的结合自由能有较大的贡献,在区分抑制剂M1(或M2)和M3的生物活性上发挥着重要的作用。     

乙酰胆碱酯酶AChE与1,7-二氮杂咔唑抑制剂的作用机理的分子动力学模拟

通过分子对接、分子动力学（MD）模拟以及成键自由能分析方法,从原子水平上模拟研究了3种1,7-二氮杂咔唑衍生物（分别记为M1、M2和M3）与AChE的结合模式及相互作用机理,分析和讨论了研究体系的静电相互作用和范德华相互作用（vdW）。用MM-PBSA方法计算的3种抑制剂与AChE之间的结合自由能与抑制剂的实验生物活性数据（IC₅₀值）相对应。分析结果表明,残基S286与抑制剂之间形成的氢键作用有利于抑制剂与AChE之间的结合。范德华相互作用,尤其是抑制剂与关键残基W279和Y334的作用,对抑制剂与AChE之间的结合自由能有较大的贡献,在区分抑制剂M1（或M2）和M3的生物活性上发挥着重要的作用。     

激酶ABL与肉豆蔻酰位点小分子作用机理的分子动力学模拟及自由能计算

采用分子动力学方法研究激酶ABL 与ATP 位点小分子imatinib、P16 及变构位点小分子STJ、MS7、MS9、3YY、MYR等的结合, 并用GBSA (generalized Born surface area)方法将结合自由能分解到各残基. 自由能计算结果表明, 小分子STJ、MS7、MS9 有利于imatinib 与ABL 结合; 小分子STJ、MS7、MS9 与激酶ABL的结合自由能接近, 绝对值均大于ABL 与3YY、MYR 的结合自由能. 能量分解表明, ABL 残基ILE502、VAL506、LEU510与STJ和MYR的相互作用是αI 螺旋处于弯曲状态的重要原因. 模拟过程中ABL肉豆蔻酰口袋残基均方根偏差(RMSD)变化值表明, STJ等小分子抑制剂与ABL结合后降低了肉豆蔻酰口袋残基的柔性.     

凝血因子Xa的S4口袋部分关键残基对利伐沙班结合的不同作用机制

通过生物信息学对比、 分子动力学模拟和结合自由能计算分析了利伐沙班与凝血因子Xa的S4口袋部分关键残基之间动态相互作用的细节. 结果表明, 利伐沙班与凝血因子Xa结合不稳定是由S4口袋关键残基突变对疏水盒子完整性的破坏所致. 其中Trp215侧链的疏水作用对抑制剂结合的作用较大, 但对整体结构的影响短时间内较小. Tyr99虽然在结合自由能中贡献较小, 但其突变可能导致99 loop所在结构域的整体构象变化, 从而对于抑制剂或底物的结合特异性产生影响. S4口袋关键残基的不同作用在凝血因子Xa直接抑制剂的药物设计及其拮抗剂的开发中应予以充分考虑.     

分子模拟方法研究四氢化吡啶并[1,2-a]吲哚酮衍生物对GSK3β和CDK5的选择性

本文通过分子对接,分子动力学模拟(MD)和MM/PBSA能量计算的方法,从分子水平研究了3个四氢化吡啶并[1,2-a]吲哚酮衍生物与CDK5和GSK3β的相互作用,并揭示了这些抑制剂对GSK3β的选择性抑制机理。分子对接结果表明,抑制剂对2种激酶具有相似的结合模式,结合口袋处的残基也都根据晶体结构的序列比对相互对应。研究体系的RMSD随时间的稳定变化,表明模拟体系已达到稳定状态,因而后续的分析是可靠的。CDK5/抑制剂体系,RMSD在0.15 nm上下波动,CDK5/M1和CDK5/M2骨架轻微波动,稍高于CDK5/M3;而GSK3β体系的RMSD值略高于CDK5体系,在0.17 nm上下波动,GSK3β/M1和GSK3ββ/M2的骨架波动平衡值则稍低于GSK3β/M3。活性较大的抑制剂增强了蛋白骨架整体的"柔性",即对激酶构象产生一定影响。能量分析表明,静电能和范德华作用能够区分不同抑制剂对同种激酶的生物活性差异。极性溶剂化自由能对区分抑制剂选择性也很重要,残基分解表明GSK3β的Glu97、Thr138是造成抑制剂选择性的主要原因。抑制剂与CDK5和GSK3β结合的过程中,蛋白质残基的动态相关性存在差异,铰链区域的Thr138与Val135~Gln206区域残基正相关,证实Thr138残基是区分抑制剂选择性的关键。     

EGFR和4-苯胺喹唑啉类抑制剂之间相互作用模式的研究

采用分子动力学和MM/PBSA相结合的方法预测了表皮生长因子受体和4-苯胺喹 啉类抑制剂的相互作用模式。在分子动力学采样的基础上,采用MM/PBSA的方法分 别预测了四种可能结合模式下表皮生长因子受体和4-苯胺喹唑啉类抑制剂间的结合 自由能。在MM/PBSA计算中,受体和抑制剂之间的非键相互作用能采用分子力学 (MM)的方法得到;溶剂效应中极性部分对自由能的贡献通过解Possion- Boltzmanne (PB)方程的方法得到;溶液效应中非极性部分对自由能的贡献则通过 分子表面积计算(SA)的方法得到。计算表明,在四种结合模式下,表皮生长因子受 体和4-苯胺喹唑啉类抑制剂之间的结合自由能有较大的差别。在最佳的相互作用模 式中,抑制剂的苯胺部分位于活性口袋的底部,能够与受体残基的非极性侧链产生 很强的范德华和疏水相互作用。抑制剂喹唑啉环上的N(1)原子能够和Met-769上的 NH形成稳定的氢键,而抑制剂上的N(3)原子则和周围的一个水分子形成氢键。同时 ,抑制剂双环上的取代基团也能和活性口袋外部的部分残基形成一定的范德华和疏 水相互作用。最佳结合模式能够很好地解释已有抑制剂结构和活性间的关系。     

FKBP12蛋白与其抑制剂结合的分子动力学模拟和结合自由能计算

FKBP12 (FK506-binding protein-12)是一种具有神经保护和促神经再生作用的蛋白. 采用分子动力学模拟取样, 运用MM-GBSA方法计算了FKBP12和3个抑制剂(GPI-1046, 308和107)的绝对结合自由能, GPI-1046的结合能最小, 308小于107的结合能. 通过能量分解的方法考察了FKBP12蛋白的主要残基与抑制剂之间的相互作用和识别, 计算结果表明: 3个抑制剂具有相似的结合模式, Ile56和Tyr82主要表现为氢键作用, Tyr26, Phe46, Val55, Ile56, Trp59, Tyr82, Tyr87和Phe99形成疏水作用区. 计算结果和实验结果吻合.     

GPCR A2A AR腺苷受体蛋白的激动剂结合及其诱导的蛋白活性开关构象变化

运用分子动力学模拟,研究了腺苷酸（激动剂）与A2AAR腺苷受体蛋白的相互作用和配体结合诱导的蛋白动力学变化．识别了与腺苷酸结合力强于0．5kcal／mol的关键基团：A63^2．61,I66^2．64,V84^3.32,L85^3.33,T88^3.36,F168^5.29,M177^5.38,L249^6.51,H250^6.52和N253^6.55,观察到腺苷酸没有与L167^5.28相互作用,这一结果支持了L167^5.28是抑制剂特异性结合位点,不与激动剂结合．未结合配体（激动剂或抑制剂）的单体A2AAR和腺苷酸结合后的A2AAR在构象上有三个不同功能性开关．腺苷酸结合可以诱导A2AAR腺苷受体蛋白的构象调整,使得三个功能性开关器件的构象与单体A2AAR不同．     

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.     

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 dynamic and free energy studies of primary resistance mutations in HIV-1 protease-ritonavir complexes

To understand the basis of drug resistance of the HIV-1 protease, molecular dynamic (MD) and free energy calculations of the wild-type and three primary resistance mutants, V82F, I84V, and V82F/I84V, of HIV-1 protease complexed with ritonavir were carried out. Analysis of the MD trajectories revealed overall structures of the protein and the hydrogen bonding of the catalytic residues to ritonavir were similar in all four complexes. Substantial differences were also found near the catalytic binding domain, of which the double mutant complex has the greatest impact on conformational changes of the protein and the inhibitor. The tip of the HIV-1 protease flap of the double mutant has the greater degree of opening with respect to that of the others. Additionally, the phenyl ring of Phe82 moves away from the binding pocket S1', and the conformational change of ritonavir subsite P1' consequently affects the cavity size of the protein and the conformational energy of the inhibitor. Calculations of binding free energy using the solvent continuum model were able to reproduce the same trend of the experimental inhibition constant. The results show that the resistance mutants require hydrophobic residues to maintain the interactions in the binding pocket. Changes of the cavity volume correlate well with free energy penalties due to the mutation and are responsible for the loss of drug susceptibility.     

Binding free energy of several sterols to the N-terminal domain of Niemann-Pick C1-like 1 protein due to mutation: Molecular dynamics study

The membrane protein Niemann-Pick type C1-like 1 (NPC1L1) plays a central role in the absorption of cholesterol in the small intestine. Other sterols, notably vitamin E and vitamin K1 also utilize NPC1L1 as a membrane transporter even though other absorption paths exist. Many NPC1L1 mutations causing the disease due to poor transport of cholesterol are known. It is not clear at this moment if the same mutation can lead to reduced transport behavior with these vitamins. In this study, we have obtained the binding free energies of these two sterols using molecular dynamics simulation and compared these values with the cholesterol-binding free energy. The N-terminal domain (NTD) of the wild as well as the disease-causing two mutations, T61M and L110F, are used for this purpose. The result indicates that the mutations show reduced binding affinity compared to the wild except for the vitamin K1 in the T61M mutant, which has increased binding free energy. Also, we found the similarity of the key amino acids responsible for the change of free energy by mutation between T61M and L110F. At the same time, non-negligible differences exist also.     

Improved Binding Free Energy Predictions from Single-Reference Thermodynamic Integration Augmented with Hamiltonian Replica Exchange

Reliable predictions of relative binding free energies are essential in drug discovery, where chemists modify promising compounds with the aim of increasing binding affinity. Conventional Thermodynamic Integration (TI) approaches can estimate corresponding changes in binding free energies, but suffer from inadequate sampling due to ruggedness of the molecular energy surfaces. Here, we present an improved TI strategy for computing relative binding free energies of congeneric ligands. This strategy employs a specific, unphysical single-reference (SR) state and Hamiltonian replica exchange (HREX) to locally enhance sampling. We then apply this strategy to compute relative binding free energies of twelve ligands in the L99A mutant of T4 Lysozyme. Besides the ligands, our approach enhances hindered rotations of the important V111, as well as V87 and L118 sidechains. Concurrently, we devise practical strategies to monitor and improve HREX-SRTI efficiency. Overall, the HREX-SRTI results agree well (R(2) = 0.76, RMSE = 0.3 kcal/mol) with available experimental data. When optimized for efficiency, the HREX-SRTI precision matches that of experimental measurements.     

Free Energy Calculations Give Insight into the Stereoselective Hydroxylation of α-Ionones by Engineered Cytochrome P450 BM3 Mutants

Previously, stereoselective hydroxylation of α-ionone by Cytochrome P450 BM3 mutants M01 A82W and M11 L437N was observed. While both mutants hydroxylate α-ionone in a regioselective manner at the C3 position, M01 A82W catalyzes formation of trans-3-OH-α-ionone products whereas M11 L437N exhibits opposite stereoselectivity, producing trans-(3S,6S)-OH-α-ionone and cis-(3S,6R)-OH-α-ionone. Here, we explore the stereoselective C3 hydroxylation of α-ionone by Cytochrome P450 BM3 mutants M01 A82W and M11 L437N using molecular dynamics-based free energy calculations to study the interaction between the enzyme and both the substrates and the products. The one-step perturbation approach is applied using an optimized reference state for substrates and products. While the free energy differences between the substrates free in solution amount to ～0 kJ mol(-1), the differences in mutant M01 A82W agree with the experimentally obtained dissociation constants K(d). Moreover, a correlation with experimentally observed trends in product formation is found in both mutants. The trans isomers show the most favorable relative binding free energy in the range of all four possible hydroxylated diastereomers for mutant M01 A82W, while the trans product from (6S)-α-ionone and the cis product from (6R)-α-ionone show highest affinity for mutant M11 L437N. Marcus theory is subsequently used to relate the thermodynamic stability to transition state energies and rates of formation.     

Converging free energy estimates: MM-PB(GB)SA studies on the protein-protein complex Ras-Raf

Estimating protein-protein interaction energies is a very challenging task for current simulation protocols. Here, absolute binding free energies are reported for the complex H-Ras/C-Raf1 using the MM-PB(GB)SA approach, testing the internal consistency and model dependence of the results. Averaging gas-phase energies (MM), solvation free energies as determined by Generalized Born models (GB/SA), and entropic contributions calculated by normal mode analysis for snapshots obtained from 10 ns explicit-solvent molecular dynamics in general results in an overestimation of the binding affinity when a solvent-accessible surface area-dependent model is used to estimate the nonpolar solvation contribution. Applying the sum of a cavity solvation free energy and explicitly modeled solute-solvent van der Waals interaction energies instead provides less negative estimates for the nonpolar solvation contribution. When the polar contribution to the solvation free energy is determined by solving the Poisson-Boltzmann equation (PB) instead, the calculated binding affinity strongly depends on the atomic radii set chosen. For three GB models investigated, different absolute deviations from PB energies were found for the unbound proteins and the complex. As an alternative to normal-mode calculations, quasiharmonic analyses have been performed to estimate entropic contributions due to changes of solute flexibility upon binding. However, such entropy estimates do not converge after 10 ns of simulation time, indicating that sampling issues may limit the applicability of this approach. Finally, binding free energies estimated from snapshots of the unbound proteins extracted from the complex trajectory result in an underestimate of binding affinity. This points to the need to exercise caution in applying the computationally cheaper "one-trajectory-alternative" to systems where there may be significant changes in flexibility and structure due to binding. The best estimate for the binding free energy of Ras-Raf obtained in this study of -8.3 kcal mol(-1) is in good agreement with the experimental result of -9.6 kcal mol(-1), however, further probing the transferability of the applied protocol that led to this result is necessary.     

Studies of Interaction Mechanism between Pyrido [3,4-d] Pyrimidine Inhibitors and Mps1

Monopolar spindle 1 (Mps1), a dual-specific kinase, is related to the proper execution of chromosome biorientation and mitotic checkpoint signaling. The overexpression of Mps1 promotes the occurrence of cancer or the survival of aneuploid cancer cells, in other words, the reduction of Mps1 will severely reduce the viability of human cancer cells. Therefore, Mps1 is a potential target for cancer treatment. Recently, a series of novel pyrido [3,4-d] pyrimidine derivatives targeting Mps1 with high biological activity were synthesized. The crystal structure of Mps1 in complex with pyrido [3,4-d] pyrimidine derivatives was also reported, but there were no specific mechanism studies for this series of small molecule inhibitors. In this study, complexes binding modes were probed by molecular docking and further validated by molecular dynamics simulations and the molecular mechanics/generalized Born surface area (MM/GBSA) method. The results indicated that the van der Waals interactions and the nonpolar solvation energies were responsible to the basis for favorable binding free energies, all inhibitors interacted with residues I531, V539, M602, C604, N606, I607, L654, I663, and P673 of Mps1. By analyzing the hydrogen bonds, we found the residues G605 and K529 in Mps1 formed stable hydrogen bonds with compounds, it was more conducive to activities of Mps1 inhibitors. According to the above analysis, we further designed five new compounds. We found that compounds IV and V were better potential Mps1 inhibitors through docking and ADMET prediction. The obtained new insights not only were helpful in understanding the binding mode of inhibitors in Mps1, but also provided important references for further rational design of Mps1 inhibitors.