Computational algorithmic and molecular dynamics study of functional and structural impacts of non-synonymous single nucleotide polymorphisms in human DHFR gene

Human dihydrofolate reductase (DHFR) is a conserved enzyme that is central to folate metabolism and is widely targeted in pathogenic diseases as well as cancers. Although studies have reported the fact that genetic mutations in DHFR leads to a rare autosomal recessive inborn error of folate metabolism and drug resistance, there is a lack of an extensive study on how the deleterious non-synonymous SNPs (nsSNPs) disrupt its phenotypic effects. In this study, we aim at discovering the structural and functional consequences of nsSNPs in DHFR by employing a combined computational approach consisting of ten recently developed in silico tools for identification of damaging nsSNPs and molecular dynamics (MD) simulation for getting deeper insights into the magnitudes of damaging effects. Our study revealed the presence of 12 most deleterious nsSNPs affecting the native phenotypic effects, with three (R71T, G118D, Y122D) identified in the co-factor and ligand binding active sites. MD simulations also suggested that these three SNPs particularly Y122D, alter the overall structural flexibility and dynamics of the native DHFR protein which can provide more understandings into the crucial roles of these mutants in influencing the loss of DHFR function.     

乙酰胆碱酯酶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的生物活性上发挥着重要的作用。     

Assessment of structurally and functionally high-risk nsSNPs impacts on human bone morphogenetic protein receptor type IA (BMPR1A) by computational approach

BMPR1A (BMP type 1 receptor) is a transmembrane cell-surface receptor also known as ALK3 (activin-like kinases-3) encodes for a type I serine/threonine kinase receptor and a member of the transforming growth-factor β–receptor (TGF-β) super family. The BMPR1A has a significant interaction with BMP-2 for protein activity and also has a low affinity with growth and differentiation factor 5 (GDF5); positively regulates chondrocyte differentiation. The genetic variations can alter the structure and function of the BMPR1A gene that causes several diseases such as juvenile polyposis syndrome or hereditary cancer-predisposing syndrome. The current study was carried out to identify potential deleterious non-synonymous single nucleotide polymorphisms (nsSNPs) in BMPR1A by implementing different computational algorithms such as SIFT, PolyPhen2, SNAP2, PROVEAN, PhD-SNP, SNPs&GO, nsSNPAnalyzer, and P-Mut. From 205 nsSNPs in BMPR1A, 7 nsSNPs (C76Y, C124R, C124Y, C376Y, R443C, R480W, and W487R) were predicted as deleterious in 8 prediction algorithms. The Consurf analysis showed that selected 7 nsSNPs were present in the highly conserved regions. Molecular dynamics simulation analysis also performed to explore conformational changes in the variant structure with respect to its native structure. According to the MDS result, all variants flexibility and rigidity were unbalanced, which may alter the structural and functional behavior of the native protein. Although, three nsSNPs i.e., C124R, C376Y, and R443C have already been reported in patients associated with JPS, but their structural and functional molecular studies remain uncharacterized. Therefore, the findings of this study can provide a better understanding of uncharacterized nsSNPS and to find their association with disease susceptibility and also facilitate to the researchers for designing or developing the target dependent drugs.     

乙酰胆碱酯酶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的生物活性上发挥着重要的作用。     

Residues R1075, D1090, R1095, and C1130 Are Critical in ADAMTS13 TSP8-Spacer Interaction Predicted by Molecular Dynamics Simulation

ADAMTS13 (A Disintegrin and Metalloprotease with Thrombospondin type 1 repeats, member 13) cleaves von Willebrand Factor (VWF) multimers to limit the prothrombotic function of VWF. The deficiency of ADAMTS13 causes a lethal thrombotic microvascular disease, thrombotic thrombocytopenic purpura (TTP). ADAMTS13 circulates in a “closed” conformation with the distal domain associating the Spacer domain to avoid off-target proteolysis or recognition by auto-antibodies. However, the interactions of the distal TSP8 domain and the Spacer domain remain elusive. Here, we constructed the TSP8-Spacer complex by a combination of homology modelling and flexible docking. Molecular dynamics simulation was applied to map the binding sites on the TSP8 or Spacer domain. The results predicted that R1075, D1090, R1095, and C1130 on the TSP8 domain were key residues that interacted with the Spacer domain. R1075 and R1095 bound exosite-4 tightly, D1090 formed multiple hydrogen bonds and salt bridges with exosite-3, and C1130 interacted with both exosite-3 and exosite-4. Specific mutations of exosite-3 (R568K/F592Y/R660K/Y661F/Y665F) or the four key residues (R1075A/D1090A/R1095A/C1130A) impaired the binding of the TSP8 domain to the Spacer domain. These results shed new light on the understanding of the auto-inhibition of ADAMTS13.     

新烟碱类化合物药效构象的分子对接研究

从烟碱型乙酰胆碱受体(nACHR)-烟碱(nicotine)复合体晶体模型出发, 采用SYBYL 6.92软件包中FlexX分子对接模块对新烟碱类化合物的3种已上市化合物吡虫啉、噻虫啉、烯啶虫胺和3种吡虫啉的结构衍生物同受体蛋白作用的精确模型进行了研究. 通过全局搜索方法构建配体的构象库进行对接, 依据构象间RMS值对结果进行分类结合CScore打分函数数据对对接结果进行筛选, 最终给出合理的新烟碱类化合物-烟碱型乙酰胆碱受体的药效作用构象模型: 配体吡啶环上氮原子通过水分子同受体Leu102, Met114形成氢键并且咪唑环或噻唑环上亲水侧链同受体CYS187或SER186形成氢键, 疏水侧链同疏水部位A (TYR164, TRP53, TYR89以及TYR185残基), 或疏水部位B (TYR132, CYS187和CYS188)相互作用. 此模型同早先有关文献报道的试验结果部分吻合, 充分表明了其合理性. 同时依据本构象模型, 在新烟碱类化合物结构方面提出了一些改良建议并为研究其高选择性指出方向.     

基于Discovery Studio软件的本科生创新实验项目设计*——卡托普利与血管紧张素转化酶的分子对接

基于Discovery Studio软件高质量的分子三维结构可视化功能和分子对接模拟模块,设计了卡托普利与血管紧张素转化酶的分子对接实验,让学生通过分子对接模拟过程,掌握药物分子结构特征,理解构效关系特征,认识酶抑制剂和靶蛋白的作用机制,了解药物研究新手段。     

基于分子对接的6-萘甲基取代HEPT类HIV-1逆转录酶抑制剂构效关系研究

采用分子对接方法得到了一系列6-萘甲基取代HEPT类逆转录酶抑制剂分子与HIV-1逆转录酶复合物模型,从中抽取出抑制剂分子的活性构象,进一步应用CoMFA和CoMSIA方法建立了具有较好预测能力的3D-QSAR模型,深入探讨了这些化合物的定量构效关系,为进一步的药物设计奠定了良好的基础.另外,以化合物13及其相应的β异构体24为代表,结合量子化学从头算分子轨道理论方法考察了它们的前线轨道,为阐明α和β系列化合物的活性差异提供了理论依据.     

雌激素β受体喹啉类配体的分子对接及3D-QSAR研究

对33个喹啉衍生物的雌激素β受体活性进行了分子对接以及比较分子力场分析(CoMFA)和比较分子相似性指数分析(CoMSIA). 对接结果显示氢键和疏水作用是配体与受体结合的主要因素,同时结果亦显示对接结合能与观测值pIC50具有极显著的线性相关性. 根据对接后各优势构象将33个样本进行叠合并进行CoMFA与CoMSIA研究,均得到了较优的结果,其中以选用立体场、静电场和疏水场建立的CoMSIA模型结果最优,其主成分数,r2,q2（LOO）和r2pred分别为2, 0.894, 0.708和0.802. 构效关系模型分析显示基团的空间位阻、电性及疏水作用是影响活性的主要因素     

Computational simulation of ligand docking to L-type pyruvate kinase subunit

Computational blind docking approach was used for mapping of possible binding sites in L-type pyruvate kinase subunit for peptides, RRASVA and the phosphorylated derivative RRAS(P_i)VA, which model the phosphorylatable N-terminal regulatory domain of the enzyme. In parallel, the same docking analysis was done for both substrates of this enzyme, phosphoenolpyruvate (PEP) and adenosine diphosphate (ADP), and for docking of fructose 1,6-bisphosphate (FBP), which is the allosteric activator of the enzyme. The binding properties of the entire surface of the protein were scanned and several possible binding sites were identified in domains A and C of the protein, while domain B revealed no docking sites for peptides or for substrates or the allosteric regulator. It was found that the docking sites of different ligands were partially overlapping, pointing to the possibility that some regulatory effects, observed in the case of L-type pyruvate kinase, may be caused by the competition of different ligands for the same binding sites.     

Molecular dynamics simulation of cocaine binding with human butyrylcholinesterase and its mutants

Molecular dynamics (MD) simulations were carried out to study cocaine binding with wild-type human butyrylcholinesterase (BChE) and its mutants based on a recently reported X-ray crystal structure of human BChE. For each BChE-cocaine system, we simulated both the nonprereactive and prereactive complexes in water. Despite the significant difference found at the acyl binding pocket, the simulated structures confirm the fundamental structural and mechanistic insights obtained from earlier computational studies of wild-type BChE with cocaine based on a homology model, e.g. the rate-determining step for BChE-catalyzed hydrolysis of biologically active (-)-cocaine is the (-)-cocaine rotation in the active site from the nonprereactive BChE-(-)-cocaine complex to the prereactive complex. It has been demonstrated that the MD simulations on both the nonprereactive and prereactive BChE-cocaine complexes can clearly reveal whether specific mutations produce the desired BChE-(-)-cocaine binding structures in which the (-)-cocaine rotation is less hindered while the required prereactive BChE-(-)-cocaine binding is maintained. Based on the MD simulations, both A328W/Y332A and A328W/Y332G BChE's are expected to have catalytic activity for (-)-cocaine hydrolysis higher than that of wild-type BChE and the activity of A328W/Y332G BChE should be slightly higher than that of A328W/Y332A BChE due to the less-hindered (-)-cocaine rotation in the mutant BChE's. However, the less-hindered (-)-cocaine rotation is only a necessary condition for a higher activity mutant BChE. The (-)-cocaine rotation is also less hindered in A328W/Y332A/Y419S BChE, but (-)-cocaine binds with A328W/Y332A/Y419S BChE in a way that is not suitable for the catalysis. Thus, A328W/Y332A/Y419S BChE is expected to lose the catalytic activity. The computational predictions were confirmed by our experimental kinetic data, demonstrating that the MD simulation-based computational protocol used in this study is reliable in prediction of the catalytic activity of BChE mutants for (-)-cocaine hydrolysis.     

In Vitro Study of the Fibrinolytic Activity via Single Chain Urokinase-Type Plasminogen Activator and Molecular Docking of FGFC1

Fungi fibrinolytic compound 1 (FGFC1) is a rare marine-derived compound that can enhance fibrinolysis both in vitro and in vivo. The fibrinolytic activity characterization of FGFC1 mediated by plasminogen (Glu-/Lys-) and a single-chain urokinase-type plasminogen activator (pro-uPA) was further evaluated. The binding sites and mode of binding between FGFC1 and plasminogen were investigated by means of a combination of in vitro experiments and molecular docking. A 2.2-fold enhancement of fibrinolytic activity was achieved at 0.096 mM FGFC1, whereas the inhibition of fibrinolytic activity occurred when the FGFC1 concentration was above 0.24 mM. The inhibition of fibrinolytic activity of FGFC1 by 6-aminohexanoic acid (EACA) and tranexamic acid (TXA) together with the docking results revealed that the lysine-binding sites (LBSs) play a crucial role in the process of FGFC1 binding to plasminogen. The action mechanism of FGFC1 binding to plasminogen was inferred, and FGFC1 was able to induce plasminogen to exhibit an open conformation by binding through the LBSs. The molecular docking results showed that docking of ligands (EACA, FGFC1) with receptors (KR1–KR5) mainly occurred through hydrophilic and hydrophobic interactions. In addition, the binding affinity values of EACA to KR1–KR5 were −5.2, −4.3, −3.7, −4.5, and −4.3 kcal/moL, respectively, and those of FGFC1 to KR1–KR5 were −7.4, −9.0, −6.3, −8.3, and −6.7 kcal/moL, respectively. The findings demonstrate that both EACA and FGFC1 bound to KR1–KR5 with moderately high affinity. This study could provide a theoretical basis for the clinical pharmacology of FGFC1 and establish a foundation for practical applications of FGFC1.     

Conjugated β-Cyclodextrin Enhances the Affinity of Folic Acid towards FRα: Molecular Dynamics Study

Drug targeting is a progressive area of research with folate receptor alpha (FRα) receiving significant attention as a biological marker in cancer drug delivery. The binding affinity of folic acid (FA) to the FRα active site provides a basis for recognition of FRα. In this study, FA was conjugated to beta-cyclodextrin (βCD) and subjected to in silico analysis (molecular docking and molecular dynamics (MD) simulation (100 ns)) to investigate the affinity and stability for the conjugated system compared to unconjugated and apo systems (ligand free). Docking studies revealed that the conjugated FA bound into the active site of FRα with a docking score (free binding energy < −15 kcal/mol), with a similar binding pose to that of unconjugated FA. Subsequent analyses from molecular dynamics (MD) simulations, root mean square deviation (RMSD), root mean square fluctuation (RMSF), and radius of gyration (Rg) demonstrated that FA and FA–βCDs created more dynamically stable systems with FRα than the apo-FRα system. All systems reached equilibrium with stable RMSD values ranging from 1.9–2.4 Å and the average residual fluctuation values of the FRα backbone atoms for all residues (except for terminal residues ARG8, THR9, THR214, and LEU215) were less than 2.1 Å with a consistent Rg value of around 16.8 Å throughout the MD simulation time (0–100 ns). The conjugation with βCD improved the stability and decreased the mobility of all the residues (except residues 149–151) compared to FA–FRα and apo-FRα systems. Further analysis of H-bonds, binding free energy (MM-PBSA), and per residue decomposition energy revealed that besides APS81, residues HIS20, TRP102, HIS135, TRP138, TRP140, and TRP171 were shown to have more favourable energy contributions in the holo systems than in the apo-FRα system, and these residues might have a direct role in increasing the stability of holo systems.     

Automated docking of ligands to an artificial active site: augmenting crystallographic analysis with computer modeling

The W191G cavity of cytochrome c peroxidase is useful as a model system for introducing small molecule oxidation in an artificially created cavity. A set of small, cyclic, organic cations was previously shown to bind in the buried, solvent-filled pocket created by the W191G mutation. We docked these ligands and a set of non-binders in the W191G cavity using AutoDock 3.0. For the ligands, we compared docking predictions with experimentally determined binding energies and X-ray crystal structure complexes. For the ligands, predicted binding energies differed from measured values by +/- 0.8 kcal/mol. For most ligands, the docking simulation clearly predicted a single binding mode that matched the crystallographic binding mode within 1.0 A RMSD. For 2 ligands, where the docking procedure yielded an ambiguous result, solutions matching the crystallographic result could be obtained by including an additional crystallographically observed water molecule in the protein model. For the remaining 2 ligands, docking indicated multiple binding modes, consistent with the original electron density, suggesting disordered binding of these ligands. Visual inspection of the atomic affinity grid maps used in docking calculations revealed two patches of high affinity for hydrogen bond donating groups. Multiple solutions are predicted as these two sites compete for polar hydrogens in the ligand during the docking simulation. Ligands could be distinguished, to some extent, from non-binders using a combination of two trends: predicted binding energy and level of clustering. In summary, AutoDock 3.0 appears to be useful in predicting key structural and energetic features of ligand binding in the W191G cavity.     

分子对接预测手性化合物在蛋白质固定相色谱柱上的保留行为

利用分子对接技术预测了2种蛋白质固定相分别与4对手性化合物的相互作用情况。 结果表明,预测结合自由能(ΔG)的大小与对映体（R-(+)型和S-(-)型）的出峰顺序一致;结合自由能差值的绝对值(Δ(ΔG))与实验分离因子(α)大小顺序一致。 说明分子对接可以反映蛋白质对不同的手性化合物的识别能力和化合物R-(+)型和S-(-)型的出峰顺序。     

Molecular Dynamics Simulations Study of the Interactions between Human Dipeptidyl-Peptidase III and Two Substrates

Human dipeptidyl-peptidase III (hDPP III) is capable of specifically cleaving dipeptides from the N-terminal of small peptides with biological activity such as angiotensin II (Ang II, DRVYIHPF), and participates in blood pressure regulation, pain modulation, and the development of cancers in human biological activities. In this study, 500 ns molecular dynamics simulations were performed on free-hDPP III (PDB code: 5E33), hDPP III-Ang II (PDB code: 5E2Q), and hDPP III-IVYPW (PDB code: 5E3C) to explore how these two peptides affect the catalytic efficiency of enzymes in terms of the binding mode and the conformational changes. Our results indicate that in the case of the hDPP III-Ang II complex, subsite S1 became small and hydrophobic, which might be propitious for the nucleophile to attack the substrate. The structures of the most stable conformations of the three systems revealed that Arg421-Lys423 could form an α-helix with the presence of Ang II, but only part of the α-helix was produced in hDPP III-IVYPW. As the hinge structure in hDPP III, the conformational changes that took place in the Arg421-Lys423 residue could lead to the changes in the shape and space of the catalytic subsites, which might allow water to function as a nucleophile to attack the substrate. Our results may provide new clues to enable the design of new inhibitors for hDPP III in the future.     

Molecular Recognition of Natural and Non-Natural Substrates by Cellodextrin Phosphorylase from Ruminiclostridium Thermocellum Investigated by NMR Spectroscopy

β-1→4-Glucan polysaccharides like cellulose, derivatives and analogues, are attracting attention due to their unique physicochemical properties, as ideal candidates for many different applications in biotechnology. Access to these polysaccharides with a high level of purity at scale is still challenging, and eco-friendly alternatives by using enzymes in vitro are highly desirable. One prominent candidate enzyme is cellodextrin phosphorylase (CDP) from Ruminiclostridium thermocellum, which is able to yield cellulose oligomers from short cellodextrins and α-d -glucose 1-phosphate (Glc-1-P) as substrates. Remarkably, its broad specificity towards donors and acceptors allows the generation of highly diverse cellulose-based structures to produce novel materials. However, to fully exploit this CDP broad specificity, a detailed understanding of the molecular recognition of substrates by this enzyme in solution is needed. Herein, we provide a detailed investigation of the molecular recognition of ligands by CDP in solution by saturation transfer difference (STD) NMR spectroscopy, tr-NOESY and protein-ligand docking. Our results, discussed in the context of previous reaction kinetics data in the literature, allow a better understanding of the structural basis of the broad binding specificity of this biotechnologically relevant enzyme.     

Molecular modelling prediction of ligand binding site flexibility

We have investigated the efficacy of generating multiple sidechain conformations using a rotamer library in order to find the experimentally observed ligand binding site conformation of a protein in the presence of a bound ligand. We made use of a recently published algorithm that performs an exhaustive conformational search using a rotamer library to enumerate all possible sidechain conformations in a binding site. This approach was applied to a dataset of proteins whose structures were determined by X-ray and NMR methods. All chosen proteins had two or more structures, generally involving different bound ligands. By taking one of these structures as a reference, we were able in most cases to successfully reproduce the experimentally determined conformations of the other structures, as well as to suggest alternative low-energy conformations of the binding site. In those few cases where this procedure failed, we observed that the bound ligand had induced a high-energy conformation of the binding site. These results suggest that for most proteins that exhibit limited backbone motion, ligands tend to bind to low energy conformations of their binding sites. Our results also reveal that it is possible in most cases to use a rotamer search-based approach to predict alternative low-energy protein binding site conformations that can be used by different ligands. This opens the possibility of incorporating alternative binding site conformations to improve the efficacy of docking and structure-based drug design algorithms.     

Assessing the performance of the molecular mechanics/Poisson Boltzmann surface area and molecular mechanics/generalized Born surface area methods. II. The accuracy of ranking poses generated from docking

In molecular docking, it is challenging to develop a scoring function that is accurate to conduct high-throughput screenings. Most scoring functions implemented in popular docking software packages were developed with many approximations for computational efficiency, which sacrifices the accuracy of prediction. With advanced technology and powerful computational hardware nowadays, it is feasible to use rigorous scoring functions, such as molecular mechanics/Poisson Boltzmann surface area (MM/PBSA) and molecular mechanics/generalized Born surface area (MM/GBSA) in molecular docking studies. Here, we systematically investigated the performance of MM/PBSA and MM/GBSA to identify the correct binding conformations and predict the binding free energies for 98 protein-ligand complexes. Comparison studies showed that MM/GBSA (69.4%) outperformed MM/PBSA (45.5%) and many popular scoring functions to identify the correct binding conformations. Moreover, we found that molecular dynamics simulations are necessary for some systems to identify the correct binding conformations. Based on our results, we proposed the guideline for MM/GBSA to predict the binding conformations. We then tested the performance of MM/GBSA and MM/PBSA to reproduce the binding free energies of the 98 protein-ligand complexes. The best prediction of MM/GBSA model with internal dielectric constant 2.0, produced a Spearman's correlation coefficient of 0.66, which is better than MM/PBSA (0.49) and almost all scoring functions used in molecular docking. In summary, MM/GBSA performs well for both binding pose predictions and binding free-energy estimations and is efficient to re-score the top-hit poses produced by other less-accurate scoring functions.