分子动力学模拟的柔性对接

分子自动对接技术在过去二十年里取得很大发展和成功,但是仍然面对如何处理分子柔性这样一个难题。这篇综述概要介绍分子柔性对接技术的进展并重点介绍分子动力学模拟技术。     

Structural transition of solvated H-Ras/GTP revealed by molecular dynamics simulation and local network entropy

The state transitions of solvated H-Ras protein with GTP were theoretically analyzed through molecular dynamics (MD) simulations. To accelerate the structural changes associated with the locations of two switch regions (I and II), the Parallel Cascade Selection MD (PaCS-MD) method was employed in this study. The interconversions between the State 1 and State 2 were thus studied in atomic details, leading to a reasonable agreement with experimental observations and consequent scenarios concerning the transition mechanism that would be essential for the development of Ras inhibitors as anti-cancer agents. Furthermore, the state-transition-based local network entropy (SNE) was calculated for the transition process from State 1 to State 2, by which the temporal evolution of information entropy associated with the dynamical behavior of hydrogen bond network composed of hydration water molecules was described. The calculated results of SNE thus proved to provide a good indicator to detect the dynamical state transition of solvated Ras protein system (and probably more general systems) from a viewpoint of nonequilibrium statistical thermodynamics.     

Effects of different solvents on the conformations of apoptotic cytochrome c: Structural insights from molecular dynamics simulation

Cytochrome c (cyt-c) upon binding with cardiolipin acquires peroxidase activity and is strictly connected to the pathogenesis of many human diseases including neurodegenerative and cardiovascular diseases. Interaction of cyt-c with cardiolipin mimics partial unfolding/conformational changes of cyt-c in different solvent environments. Dynamic pictures of these conformational changes of cyt-c are crucial in understanding their physiological roles in mitochondrial functions. Therefore, atomistic molecular dynamics (MD) simulations have been carried out to investigate the effect of different solvents (water, urea/water, MeOH and DMSO) on the structure and conformations of apoptotic cyt-c (Fe³⁺). Our study demonstrates that the structural changes in the protein are solvent dependent. The structural differences are observed majorly on the β-sheets and α-helical conformations and the degree of their perturbation are specific to the solvent. Although a complete loss of β-sheets (0%) is observed in MeOH and DMSO, by contrast, well preserved β-sheets (3.84%) are observed in water and urea/water. A significant decrease in the α-helical contents is observed in MeOH (41.34%) and water (42.46%), a negligible alteration in DMSO (44.25%) and well preserved α-helical (45.19%) contents in urea/water. The distances between the residues critical for electron transfer are decreased considerably for DMSO. Further, the reduction in residue flexibility and the conformational space indicate that the collective motions of cyt-c are reduced when compared to other cosolvents. Essential dynamics analysis implies that the overall motions of cyt-c in water, MeOH and urea/water are involved in three to four eigenvectors and in first eigenvector in DMSO. Overall, we believe that MD simulations of cyt-c in different solvents can provide a detailed microscopic understanding of the physiological roles, electron transport and peroxidase function in the early events of apoptosis which are hard to probe experiments.     

Algorithm optimization in molecular dynamics simulation

Establishing the neighbor list to efficiently calculate the inter-atomic forces consumes the majority of computation time in molecular dynamics (MD) simulation. Several algorithms have been proposed to improve the computation efficiency for short-range interaction in recent years, although an optimized numerical algorithm has not been provided. Based on a rigorous definition of Verlet radius with respect to temperature and list-updating interval in MD simulation, this paper has successfully developed an estimation formula of the computation time for each MD algorithm calculation so as to find an optimized performance for each algorithm. With the formula proposed here, the best algorithm can be chosen based on different total number of atoms, system average density and system average temperature for the MD simulation. It has been shown that the Verlet Cell-linked List (VCL) algorithm is better than other algorithms for a system with a large number of atoms. Furthermore, a generalized VCL algorithm optimized with a list-updating interval and cell-dividing number is analyzed and has been verified to reduce the computation time by 30∼60% in a MD simulation for a two-dimensional lattice system. Due to similarity, the analysis in this study can be extended to other many-particle systems.     

Scalability study of molecular dynamics simulation on Godson-T many-core architecture

Molecular dynamics (MD) simulation has broad applications, and an increasing amount of computing power is needed to satisfy the large scale of the real world simulation. The advent of the many-core paradigm brings unprecedented computing power, but it remains a great challenge to harvest the computing power due to MD’s irregular memory-access pattern. To address this challenge, this paper presents a joint application/architecture study to enhance the scalability of MD on Godson-T-like many-core architecture. First, a preprocessing approach leveraging an adaptive divide-and-conquer framework is designed to exploit locality through memory hierarchy with software controlled memory. Then three incremental optimization strategies–a novel data-layout to improve data locality, an on-chip locality-aware parallel algorithm to enhance data reuse, and a pipelining algorithm to hide latency to shared memory–are proposed to enhance on-chip parallelism for Godson-T many-core processor. Experiments on Godson-T simulator exhibit strong-scaling parallel efficiency of 0.99 on 64 cores, which is confirmed by a field-programmable gate array emulator. Also the performance per watt of MD on Godson-T is much higher than MD on a 16-cores Intel core i7 symmetric multiprocessor (SMP) and 26 times higher than MD on an 8-core 64-thread Sun T2 processor. Detailed analysis shows that optimizations utilizing architectural features to maximize data locality and to enhance data reuse benefit scalability most. Furthermore, a hierarchical parallelization scheme is designed to map the MD algorithm to Godson-T many-core cluster and a simple performance model is derived, which suggests that the optimization scheme is likely to scale well toward exascale. Certain architectural features are found essential for these optimizations, which could guide future hardware developments.     

Huge-scale molecular dynamics simulation of multibubble nuclei

We have developed molecular dynamics codes for a short-range interaction potential that adopt both the flat-MPI and MPI/OpenMP hybrid parallelizations on the basis of a full domain decomposition strategy. Benchmark simulations involving up to 38.4 billion Lennard-Jones particles were performed on Fujitsu PRIMEHPC FX10, consisting of 4800 SPARC64 IXfx 1.848 GHz processors, at the Information Technology Center of the University of Tokyo, and a performance of 193 teraflops was achieved, which corresponds to a 17.0% execution efficiency. Cavitation processes were also simulated on PRIMEHPC FX10 and SGI Altix ICE 8400EX at the Institute of Solid State Physics of the University of Tokyo, which involved 1.45 billion and 22.9 million particles, respectively. Ostwald-like ripening was observed after the multibubble nuclei. Our results demonstrate that direct simulations of multiscale phenomena involving phase transitions from the atomic scale are possible and that the molecular dynamics method is a promising method that can be applied to petascale computers.     

Structural insights into the camel milk lactoperoxidase: Homology modeling and molecular dynamics simulation studies

Lactoperoxodase (LPO) is a heme peroxidase enzyme present in mammalian milk. It is an antimicrobial protein with wide range of industrial applications. Although the three dimensional structure of LPO from various mammalian species has been reported, but its structure from camel source is still unknown. So far, the crystallization attempts have not been successful in determining camel LPO (cLPO) structure. Herein, we developed the three dimensional structure of cLPO by homology modeling approach using prime module available in Schrodinger suite. The developed model in complex with ligand hypothiocyanate (OSCN⁻) was further validated by Ramachandran plot followed by molecular dynamics (MD) simulation studies using Desmond module of Schrodinger. cLPO model exhibited overall structural similarity with template crystal structure, however, it displayed different interaction pattern of amino acid residues with ligand OSCN⁻ in comparison to template crystal structure. Moreover, the ligand binding site environment in cLPO is more polar, less hydrophobic, and harbours more number of charged residues than template crystal structure. The substrate binding pocket environment of cLPO shows a considerable difference from template crystal structure. This subsequently resulted in dissimilar behaviour of ligand during the course of MD simulation studies.     

微孔隙中流体扩散系数分子动力学模拟的并行算法研究

扩散系数计算量非常大,为提高计算效率研究采用并行计算替代单机计算：用非阻塞通信替代阻塞通信;利用元胞分解的成果创建要拷贝的外层粒子列表;用拷贝粒子和移动粒子的重叠部分来减少通信量,在空间分解(SD)算法的基础上对MD模拟的并行算法进行了改进,使并行计算的效率提高了28．64％。用改进的SD并行计算方法,以RedHat Linux 8．0位操作系统和基于Message Passing Interface(MPI)消息传递界面的四节点计算机集群系统(PC Cluster)。对水分子在宏量条件下以及受限于硅酸岩狭缝中的速度自相关函数进行了MD模拟计算。结果表明,并行算法可以使速度自相关函数及相应的自扩散系数的MD模拟时间大为缩短。     

Insights into thermal stability of thermophilic nitrile hydratases by molecular dynamics simulation

Thermal stability is of great importance for industrial enzymes. Here we explored the thermal-stable mechanism of thermophilic nitrile hydratases (NHases) utilizing a molecular dynamic simulation. At a nanosecond timescale, profiles of root mean square fluctuation (RMSF) of two thermophilic NHases, 1UGQ and 1V29, under enhancing thermal stress were carried out at 300 K, 320 K, 350 K and 370 K, respectively. Results showed that the region A1 (211-231 aa) and A2 (305-316 aa) in 1UGQ, region B1 (186-192 aa) in 1V29, and most of terminal ends in both enzymes are hyper-sensitive. Salt-bridge analyses revealed that in one hand, salt-bridges contributed to maintaining the rigid structure and stable performance of the thermophilic 1UGQ and 1V29; in the other hand, salt-bridges involved in thermal sensitive regions are relatively weak and prone to be broken at elevated temperature, thereby cannot hold the stable conformation of the spatial neighborhood. In 1V29, region A1 was stabilized by a well-organized hook-hook like cluster with multiple salt-bridge interactions, region A2 was stabilized by two strong salt-bridge interactions of GLU52-ARG332 and GLU334-ARG332. In 1UGQ, the absence of a charged residue decreased its thermal sensitivity of region B1, and the formation of a small beta-sheet containing a stable salt-bridge in C-beta-terminal significantly enhanced its thermal stability. By radius of gyration calculation containing or eliminating the thermal sensitive regions, we quantified the contribution of thermal sensitive regions for thermal sensitivity of 1UGQ and 1V29. Consequently, we presented strategies to improve thermal stability of the industrialized mesophilic NHase by introducing stable salt-bridge interactions into its thermal sensitive regions.     

分子动力学并行算法的优化与应用

比较了小型机群中模拟小规模粒子系统采用原子分解和区域分解算法的模拟时间.针对原子分解算法,优化了各进程粒子的受力计算,尤其对严重影响并行效率的全局通信进行了优化,最终使得并行效率提高了30%以上.使用该优化算法对固态氩的导热系数进行了计算,模拟结果与实验值比较吻合,大大缩短了时间.     

Computational study of human phosphomannose isomerase: Insights from homology modeling and molecular dynamics simulation of enzyme bound substrate

Phosphomannose isomerase is a zinc metalloenzyme that catalyzes the reversible isomerization of mannose-6-phosphate and fructose-6-phosphate, and the three-dimensional (3D) structure of human phosphomannose isomerase has not been reported. In order to understand the catalytic mechanism, the 3D structure of the protein is built by using homology modeling based on the known crystal structure of mannose-6-phosphate isomerase from (PDB code 1PMI). The model structure is further refined by energy minimization and molecular dynamics methods. The mannose-6-phosphate-enzyme complex is developed by molecular docking and the key residues involved in the ligand binding are determined, which will facilitate the understanding of the action mode of the ligands and guide further genetic studies. Our results suggest a hydride transfer mechanism of alpha-hydrogen between the C1 and C2 positions but do not support the cis-enediol mechanism. The detailed mechanism involves, on one side, Zn2+ mediating the movement of a proton between O1 and O2, and, on the other side, the hydrophobic environment formed in part by Tyr278 promoting transfer of a hydride ion.     

Membrane protein structure quality in molecular dynamics simulation

Our goal was to assess the relationship between membrane protein quality, output from protein quality checkers and output from molecular dynamics (MD) simulations. Membrane transport proteins are essential for a wide range of cellular processes. Structural features of integral membrane proteins are still under-explored due to experimental limitations in structure determination. Computational techniques can be used to exploit biochemical and medium resolution structural data, as well as sequence homology to known structures, and enable us to explore the structure-function relationships in several transmembrane proteins. The quality of the models produced is vitally important to obtain reliable predictions. An examination of the relationship between model stability in molecular dynamics (MD) simulations derived from RMSD (root mean squared deviation) and structure quality assessment from various protein quality checkers was undertaken. The results were compared to membrane protein structures, solved at various resolution, by either X-ray or electron diffraction techniques. The checking programs could predict the potential success of MD in making functional conclusions. MD stability was shown to be a good indicator for the quality of structures. The quality was also shown to be dependent on the resolution at which the structures were determined.     

A general purpose parallel molecular dynamics simulation program

We present a general purpose parallel molecular dynamics simulation code. The code can handle NVE, NVT, and NPT ensemble molecular dynamics, Langevin dynamics, and dissipative particle dynamics. Long-range interactions are handled by using the smooth particle mesh Ewald method. The implicit solvent model using solvent-accessible surface area was also implemented. Benchmark results using molecular dynamics, Langevin dynamics, and dissipative particle dynamics are given.

Program summary

Title of program:MM_PARCatalogue identifier:ADXP_v1_0Program summary URL:http://cpc.cs.qub.ac.uk/summaries/ADXP_v1_0Program obtainable from: CPC Program Library, Queen's University of Belfast, N. IrelandComputer for which the program is designed and others on which it has been tested:any UNIX machine. The code has been tested on Linux cluster and IBM p690Operating systems or monitors under which the program has been tested:Linux, AIXProgramming language used:CMemory required to execute with typical data:∼60 MB for a system of atoms Has the code been vectorized or parallelized? parallelized with MPI using atom decomposition and domain decompositionNo. of lines in distributed program, including test data, etc.:171 427No. of bytes in distributed program, including test data, etc.:4 558 773Distribution format:tar.gzExternal routines/libraries used:FFTW free software (http://www.fftw.org)Nature of physical problem:Structural, thermodynamic, and dynamical properties of fluids and solids from microscopic scales to mesoscopic scales.Method of solution:Molecular dynamics simulation in NVE, NVT, and NPT ensemble, Langevin dynamics simulation, dissipative particle dynamics simulation.Typical running time:Table below shows the typical run times for the four test programs.

Benchmark results. The values in the parenthesis are the number of processors used
System	Method	Timing for 100 steps in seconds
256 TIP3P	MD	23.8 (1)
64 DMPC + 1645 TIP3P	MD	890 (1)	528 (2)	326 (4)	209 (8)
8 Aβ16-22	LD	1.02 (1)
23760 Groot-Warren particles	DPD	22.16 (1)

Full-size table

     

Comparative structural studies of psychrophilic and mesophilic protein homologues by molecular dynamics simulation

Comparative molecular dynamics simulations of psychrophilic type III antifreeze protein from the North-Atlantic ocean-pout Macrozoarces americanus and its corresponding mesophilic counterpart, the antifreeze-like domain of human sialic acid synthase, have been performed for 10 ns each at five different temperatures. Analyses of trajectories in terms of secondary structure content, solvent accessibility, intramolecular hydrogen bonds and protein–solvent interactions indicate distinct differences in these two proteins. The two proteins also follow dissimilar unfolding pathways. The overall flexibility calculated by the trace of the diagonalized covariance matrix displays similar flexibility of both the proteins near their growth temperatures. However at higher temperatures psychrophilic protein shows increased overall flexibility than its mesophilic counterpart. Principal component analysis also indicates that the essential subspaces explored by the simulations of two proteins at different temperatures are non-overlapping and they show significantly different directions of motion. However, there are significant overlaps within the trajectories and similar directions of motion of each protein especially at 298 K, 310 K and 373 K. Overall, the psychrophilic protein leads to increased conformational sampling of the phase space than its mesophilic counterpart.Our study may help in elucidating the molecular basis of thermostability of homologous proteins from two organisms living at different temperature conditions. Such an understanding is required for designing efficient proteins with characteristics for a particular application at desired working temperatures.     

丁二酰亚胺类分散剂在伪烟炱表面吸附行为的分子动力学模拟

丁二酰亚胺类分散剂是一种典型的用于各种润滑油中的无灰分散剂。由于烟炱颗粒结构的复杂性,从实验研究的角度对分散剂的作用机理尚未给出满意的解释。本文应用分子动力学模拟方法,研究了丁二酰亚胺类分散剂模型分子在伪烟炱表面的吸附行为。模拟结果可以为分散剂的分子设计提供参考。     

Crystal morphology prediction of 1,3,3-trinitroazetidine in ethanol solvent by molecular dynamics simulation

In order to understand the mechanism of the effect of solvent on the crystal morphology of explosives, and be convenient for the choice of crystallization solvent, the attachment energy (AE) model was performed to predict the growth morphology and the main crystal faces of 1,3,3-trinitroazetidine (TNAZ) in vacuum. The molecular dynamics simulation was applied to investigate the interactions of TNAZ crystal faces and ethanol solvent, and the growth habit of TNAZ in ethanol solvent was predicted using the modified AE model. The results indicate that the morphology of TNAZ crystal in vacuum is dominated by the six faces of [0 2 1], [1 1 2], [0 0 2], [1 0 2], [1 1 1] and [0 2 0], and the crystal shape is similar to polyhedron. In ethanol solvent, The binding strength of ethanol with TNAZ faces changes in the order of [0 2 1] > [1 1 2] > [0 0 2] > [1 0 2] > [1 1 1] > [0 2 0], which causes that [1 1 1] and [0 2 0] faces disappear and the crystal morphology becomes more regular. The radial distribution function analysis shows that the interactions between solvent and crystal faces mainly consist of coulomb interaction, van der Waals force and hydrogen bonds.     

Comparative structural studies of psychrophilic and mesophilic protein homologues by molecular dynamics simulation

Comparative molecular dynamics simulations of psychrophilic type III antifreeze protein from the North-Atlantic ocean-pout Macrozoarces americanus and its corresponding mesophilic counterpart, the antifreeze-like domain of human sialic acid synthase, have been performed for 10 ns each at five different temperatures. Analyses of trajectories in terms of secondary structure content, solvent accessibility, intramolecular hydrogen bonds and protein–solvent interactions indicate distinct differences in these two proteins. The two proteins also follow dissimilar unfolding pathways. The overall flexibility calculated by the trace of the diagonalized covariance matrix displays similar flexibility of both the proteins near their growth temperatures. However at higher temperatures psychrophilic protein shows increased overall flexibility than its mesophilic counterpart. Principal component analysis also indicates that the essential subspaces explored by the simulations of two proteins at different temperatures are non-overlapping and they show significantly different directions of motion. However, there are significant overlaps within the trajectories and similar directions of motion of each protein especially at 298 K, 310 K and 373 K. Overall, the psychrophilic protein leads to increased conformational sampling of the phase space than its mesophilic counterpart.Our study may help in elucidating the molecular basis of thermostability of homologous proteins from two organisms living at different temperature conditions. Such an understanding is required for designing efficient proteins with characteristics for a particular application at desired working temperatures.     

多核并行技术在分子动力学模拟中的应用

为了充分利用多核处理器资源,研究了一种用于分子动力学模拟中的多核并行技术。在多核处理器上利用OpenMP技术实现多线程创建与同步、动态设置子线程的调度运行方式以及负载均衡以减少子线程执行等待时间。通过对不同分子体系结构下的动力学模型测试,得出在不同子线程下并行计算的时间,并且得到了良好的性能加速比。实验结果表明,采用OpenMP并行技术可有效地提高电荷求解过程在分子动力学模拟运算中的时间效率,以及多核计算机资源的利用率。     

烷基芳基磺酸盐界面性质的分子动力学模拟

利用分子动力学方法模拟了阴离子表面活性剂在水/辛烷,水/甲苯,水/辛烷、甲苯混合油相时的动力学特点,结果很好地展示了体系静止时清楚的油水界面,并与无表面活性剂的油水界面进行了比较.利用密度剖面图可说明表面活性剂对不同油相的诱导作用.同时,模拟计算出烷基芳基磺酸盐表面活性剂在油水界面的界面张力、均方根偏差、接触面积等,由计算结果可知,油相为甲苯时能很好的降低其界面张力,油相为辛烷时更易于形成清楚的油水界面且界面强度较大.     

Parallel implementation of molecular dynamics simulation for short-ranged interaction

A parallel molecular dynamics simulation method, designed for large-scale problems, employing dynamic spatial domain decomposition for short-ranged molecular interactions is proposed. In this parallel cellular molecular dynamics (PCMD) simulation method, the link-cell data structure is used to reduce the searching time required for forming the cut-off neighbor list as well as for domain decomposition, which utilizes the multi-level graph-partitioning technique. A simple threshold scheme (STS), in which workload imbalance is monitored and compared with some threshold value during the runtime, is proposed to decide the proper time for repartitioning the domain. The simulation code is implemented and tested on the memory-distributed parallel machine, e.g., PC-cluster system. Parallel performance is studied using approximately one million L-J atoms in the condensed, vaporized and supercritical states. Results show that fairly good parallel efficiency at 49 processors can be obtained for the condensed and supercritical states (∼60%), while it is comparably lower for the vaporized state (∼40%).