gromacs文件介绍and一些杂知识

gromacs(GMX) 各种文件格式详细，可以查阅GROMACS 手册第5章第6小节，以下为简要介绍。 （1）gromacs(GMX) 各种文件格式详细，可以查阅GROMACS 手册第5章第6小节，以下为简要介绍。 CPT文件：该文件为模拟断点文件（check point，.cpt）。该文件为模拟过程固定时间间隔产生，保存模拟系统所有信息。该文件一部分可以在能量文件（.edr）找到，一部分可以在双精度轨迹文件（.trr）中找到。如果模拟不幸因为外界条件中断（如断电，模拟人发脾气砸电脑等），可以使用该文件重新在断点处开始模拟，以节省模拟时间。同时也可以依靠该断点文件开始，并延长模拟计算（见tpbconv）。 EDR文件：系统能量文件（energy，.edr）。该文件记录模拟输入文件中定义的能量组的各种相互作用能量等。 EPS文件：封装文件格式（.eps），并不是GROMACS自身文件格式，可以当图片打开。LINUX系统下一般已经有默认打开程序，WINDOWS要安装其他打开程序（可以GOOGLE以下）。GROMACS的DSSP和罗麽占陀罗图等通过xpm2ps处理后都是这个文件格式。习惯就好。 G87文件：分子坐标文件（.g87）。该文件记录并只记录原子坐标和速度，不含原子序号。并只记录常压强模拟系统的盒子信息。 G96文件：分子坐标文件（.g96）。GROMOS96程序的分子坐标文件，模拟程序以15.9的C语言格式写入，精度较高，但是会比较大。包含有文件头，时间步，原子坐标，原子速度，以及盒子信息等。 GRO文件：分子坐标文件（.gro）。GROMACS的最主要分子坐标文件，明白这个文件，就基本明白使用GROMACS了。该文件类型的各个文本列字数固定，C语言的写入格式为："%5d%5s%5s%5d%8.3f%8.3f%8.3f%8.4f%8.4f%8.4f"。具体固定文本列有：残基序号，5位数；残基名称，5字母；原子名称，5字母；原子序号，5为数；原子坐标三列，X，Y，Z坐标各8位数，含3个小数位；速度同坐标，速度单位为nm/ps（km/s）。 ITP文件：分子拓扑文件（.itp）。被主拓扑文件（.top）包含的分拓扑文件，一般包含某个特定分子的类型。于主拓扑文件区别有它不引用其他力场文件，同时包含[system]，[molecule]等拓扑字节。 M2P文件：xpm2ps程序配置文件，定义输出eps文件中颜色，字体种类及大小等。 MDP文件：GROMACS的模拟配置文件（.mdp）。该文件所含定义较多，各关键字的含义可以查阅GROMACS手册。（这是使用GROMACS进行分子动力学模拟最最最最（10个最）重要的文件，no mdp文件，no GROMACS模拟。好好看书，以明白各个关键字的含义。因为它太重要，所以不在此简要描述。 N2T文件：原子名称及类型对照文件（.n2t）。x2top程序可以按照原子名称得到该原子的原子类型力场参数，N2T就是x2top程序扫描的数据库，文件很小。文件中文本行有原子名称，原子类型，原子电量，原子质量，该原子与其他原子成键距离等。 NDX文件：原子索引文件（.ndx）。该文件含原子的序号，当使用make_ndx程序生成索引文件时，可以定义不同的原子组，每组名下即是该组所含各个原子的序号。 PDB文件：分子坐标文件（.pdb）。这个就不用说了（说真的，如果真没有听过这个文件类型的话，看这篇文章有点浪费时间。） RTP文件：残基力场参数文件（.rtp）。该文件包含常见残基的力场信息，包括残基所含原子，成键种类等。使用pdb2gmx处理PDB文件时，程序按照PDB文件信息，在RTP文件中寻找对应的残基力场信息。 TOP文件：模拟系统的拓扑文件（.top）。该文件就是所谓十分及其著名的系统拓扑文件啦，其包含各个关键字都十分易懂；一般其还包含引用其他力场文件（#include）。TOP文件一般由pdb2gmx产生，grompp程序生成模拟TPR文件时使用。 TPR文件：模拟打包文件（.tpr）。该文件打包模拟需要各种信息，包括模拟系统，模拟控制等。 TRJ文件：全精度轨迹文件（.trj）。该文件包含模拟系统模拟各个时间下的原子坐标，速度和受力等。所含帧数频率由MDP文件控制，文件较大。 TRR文件：以上同，一般为默认格式。由于所含信息多，可以也EDR文件一起使用，重新开始模拟程序。 XPM文件：数据矩阵文件（.xpm）。该文件矩阵中每个值即是矩阵点所表示的物理量大小（也可以是布尔值）。该文件其实就是二维图，可以失踪xpm2ps转换为图片。 XTC文件：模拟轨迹单精度文件（.xtc）。单精度轨迹文件，文件较TRR和TRJ小，为常用分析文件。包含模拟系统中原子坐标，模拟时间，和模拟盒子信息。 XVG文件：二维图标文件（.xvg）。二维画图工具xmgrace的默认文件，可以使用xmgrace打开。 （2）Gromacs中几个特殊文件 aminoacids.dat 该文件保存GMX默认的蛋白质和核算的默认残基名称。如果计算过程要建立一个新的蛋白质或者核算残基，可以将新的残基名称加到该文件中，并增加文件第一个的整数即可。有时候可以将该文件拷贝到当前工作文件夹进行编辑，以不影响其他计算的命名（GMX的文件搜索总是从当前目录开始的。） FF.dat GMX默认力场列表，即pdb2gmx处理PDB文件时可以选择的立场列表。增加新的力场，可以编辑该文件，并修改文件第一行的整数，使其与力场种类熟目一致。 specbond.dat GMX处理特殊化学键的文件，特殊化学键包括二硫键，血红素铁原子于其他原子成键等。该文件第一行指明特殊键对的数目，第二行开始即为各个特殊键对的信息，其中第一列为键对第一个残基的名称，第二列为该残基成键原子的名称，第三列为该原子可以成键的数目，第四到第六列为成键另一个残基的信息，第七列为该化学键的平衡长度，此后两列为成键后残基的新名称。 vdwradii.dat 原子范德华半径数据库。使用genbox为系统添加水分子，或者使用genion为系统添加离子时，各个原子间的距离要大于两个原子范德华半径之和，否则则为原子重叠 （3）常见水分子模型 进行分子动力学模拟，水分子十分重要，除非选择使用连续介质模型（implictit water model）。水分子模型较多，选择这些模型要结合使用的力场，并参考别人已经的数据。一下简单介绍几种常见的水分子模型，希望对了解它们有点帮助。 按照一般化学常识，水分子由三个原子构成，主要的参数应该有各个原子的质量，电量，氢氧键的长度以及H-O-H的键角。没有错，最简单的水分子模型就是这些参数都固定的刚性水分子模型。如SPC模型和TIP3P模型。这两种模型中，原子质量和电量都在同一个质点上。唯一不同的是TIP3P的H-O-H键角比理论值109.47小，为104.52度。这两种水模型只有氧原子具有范德华作用系数，氢原子的范德华系数为0。 以上两种模型有对应的改进模型，SPC的改进模型为SPC/E，起主要改进其实就是使溶液系统的总能量乘以5.22 kJ/mol。这样可以使SPC溶液属性更加接近实验值。TIP3P在CHARMM力场中的改进是给氢原子一定的范德华系数，这样做的结果的计算根据复杂。。。（很无奈，因为结果好，所以也没有 办法 鲁班奖评选办法下载鲁班奖评选办法下载鲁班奖评选办法下载企业年金办法下载企业年金办法下载 。） 由于真是情况下水分子的电量分布并不是完全在原子上的，如氧原子的一部分负电量就在H-O-H的对角线上，还有两个电子对处在H-O化学键的延长线上。为了得到更加真实的水分子模型，四个粒子以上的模型就被应用到分子动力学模拟中。其中最著名的有TIP4P模型。该模型在三个原子中间，H-O-H化学键的对角线上多了一个不含质量，只带电量的点。很多蛋白质模拟计算中，TIP4P和OPLS力场结合使用都得到很好的效果。 以上提到，水分子的氧原子在H-O化学键延长线上有两个电子对，于是有的人就在这两处添加了两个只带电量的粒子。2000年报道的TIP5P模型，计算结果也很好。还有一些牛人，结合TIP4P和TIP5P，要研制TIP6P，很好很强大。。。 不得不说，并不是模型的所含粒子越多越好。粒子越多，就算付出越大，因为要计算的相互作用更多 （4）力场 “力场”，请不要被“场”这个听起来像是十分高深的物理名词给吓坏了。 分子动力学模拟中使用的力场，包含两个重要的部分： 1）模拟粒子之间相互作用的方程（即经典力学的相互作用力方程如库仑定律，范德华作用方程等）。 2）方程的参数（即各个不同粒子，原子本身的参数，如带点量等等）。 可以想想，模计算机模拟好多成键或者不成键的粒子的运动，总要让它们互相推推拉拉吧，于是力场就是定义它们推推拉拉的方式（按照物理定律）。 力场类型，一般分类为三种： i）全原子力场：精确定义每一个原子的参数。 ii）联合原子力场：省略非极性氢原子，同时把其参数整合到与他们成键的相邻原子上（比如甲基，只由一个碳原子表示）。 iii）粗颗粒力场：进一步精简分子结构的力场参数，种类比较多，比如有讲蛋白侧链看作一个颗粒的力场，或者甚至将整个氨基酸残基看成一个颗粒的力场等等。 一般来说力场的方程和参数是自成一个系统的，所以一般不能在一个系统中使用两个力场的参数。更具体的将，同一个原子在力场一中的带电量与起在第二个力场中是不一样的，化学键也一样。一般来讲，也不能特定修改力场中模一个原子的参数，因为原子之间是互相交叠依赖（比如未来保证整个氨基酸残基电量为0，各个原子电量加和必须为0）。但是，这并不是说一定不行，相反的，为了模拟一些不常见的分子，经常需要根据已有的参数（力场里面的，其他论文等）来构建新的分子参数。具体方法可以参考Mr. Google等著名老师。 目前比较流行的力场有： AMBER： 包含好几个版本的力场，为全原子力场； CHARMM：全原子力场，是软件CHARMM的一部分； GROMOS：GROMOS软件使用的力场，版本较多，为联合原子力场；OPLS：包含全原子和联合原子力场两个版本； 粗颗粒力场：种类较多，没有固定版本或者种类，一般根据研究需要开发。 （5）Gromacs重启模拟计算 以前介绍过如果使用 GMX 3.x 重新由于种种原因停止的模拟，以下为 GMX 4.x 下重启模拟的方法。 GMX 4.x 的模拟程序 mdrun 较以往版本有不少不同。在模拟过程中，mdrun 按照 mdp 文件在一定时间间隔保存一个断点文件（checkpoint file, .cpt文件），该文件保存了该时刻模拟系统的所有物理量信息。如果由于不可预见原因，模拟中断，则可以使用该文件重新在该时刻开始进行模拟。 重启模拟的命令如下： ------- mdrun -s topol.tpr -cpi state.cpt -append ------- 以上 state.cpt文件为最新生产的断点文件（ mdrun 会保存另外一个断点文件：state_prev.cpt，为上一个时刻保存的断点文件，双保险。）使用 “-append " 的作用是将模拟输出添加到已有文件中，包括轨迹文件，记录文件，能量文件等，相同帧的信息将被后生产的信息覆盖。 当然，也可以继续像 GMX 3.x 一样使用 tpbconv生产新的 tpr 文件继续模拟，详细请参见旧文或手册。 （6）Gromacs多链模拟 进行模拟计算时，如果模拟分子由两条以上的链组成，一般都要明确告诉模拟软件区分两条链。模拟软件一般没有那么聪明，除非明确定义，否则它会把两条以上的化学链（如肽链，DNA，其他聚酰胺等）看成一条链。在建立模拟文件是，上一条链尾端会于下一条链头部加一个共价化学键（如肽键）。由于该化学键一般很长，开始模拟时系统就“爆炸”了。 AMBER软件在处理这样的问题的，需要编辑原始的PDB文件，在每一条链结尾处添加“TER”。在GMX中，这种做法行不通（其实开发人员应该考虑这个问题）。解决的办法要在原始PDB文件中给每一条链添加链标识符，如“A”，“B”等等。（如果26个字母不够用，那就使用数字1到9，然后还可以使用特殊字符，如"$“，”￥“等等）。这样，使用 pdb2gmx 处理PDB文件的时候，就会得到各个链的拓扑文件，如topol_A.itp，topol_B.itp等等，并都被topol.top包含。 以上所述使用一个字符标识PDB文件中不同的链，是因为PDB文件只使用第22字符列作为链标识位，两个字符以上不认。（即AA，AB标识的链都被认为是A链。） 那么如果拿到一个没有链标识符的PDB坐标文件或GRO文件，该怎么办呢？那么要先使用make_ndx将不同链的残基选作不同的分子组（group），然后使用editconf将不同组输出成带链标识符的PDB文件，命令如： editconf -f File.pdb(/File.gro) -n indenx.ndx -o chian_A.pdb -lable A （以上可以等等A链，以此类推得到不同的链的PDB文件）。 最后将这些PDB文件组合成一个PDB文件，再由pdb2gmx处理即可。甚是麻烦。。。 （7）GMX.5 eneconv GMX 分子模拟，有一个非常重要的能量输出文件，即 edr 文件。eneconv 就是对 GMX 能量输出文件进行处理的程序。 一 个模拟可以分对次进行，于是得到很多 edr 文件。使用 eneconv 的 “ -f ” 参数，然后把这些能量文件罗列出来，那么就可以对这些能量文件进行合并，输出一个完整的能量文件。如果另个能量文件中有重复的模拟步骤，那么后一个读入的 能量文件将覆盖前一个文件。也可以使用 “-settime” 参数对每一个输入文件的开始时间进行设置，以免互相覆盖。如下就是一个程序运行例子： ++++++++++++++++ eneconv -o fixed.edr -f *.edr ++++++++++++++++ 即对当前目录下所有 edr 文件进行合并，然后输出为 fixed.edr文件。 当用 “ -f ” 参数读入单一一个能量 edr 文件时，也不是没有用，可以和其他参数配合，对能量文件进行编辑, 如 “ -dt ” 参数可以设定对原来能量文件进行规定时间间隔输出到新能量文件中；“ -offset ” 参数设定写出输入能量文件的时间帧（从哪一个模拟时间开始写入新的能量文件）等。这些参数，还有上面说到的 “-settime” 参数，都可以一起使用，加上 “-b” 和 “-e” 设定开始和结束读取模拟时间帧，就能得到新的称心如意的新能量文件。 程序输入文件： ----------- -f: 输入能量文件，即 edr 文件。 -o: 输出文件，也是 edr 文件。 其他参数： --------- -b: 设定从哪一个模拟时间帧对输入文件进行读取。 -e: 设定从哪一个模拟时间帧对输入文件结束读取。 -dt: 设定输出文件的模拟时间间隔，比如 “-dt 10” 表示每 10ps 输出一次。 -offset: 设定从哪一个时间帧开始输出到新的能量文件中。 -settime: 交换式设定每一个输入文件在新输出文件中的开始时间。 -sort: 自动排序输入文件。 -scalefac: 该参数输入为一个实数，程序会将能量文件中的每一个能量项乘以这一个实数。 -error: 如果输入文件中有错误，程序自动退出。 注意：新的输出文件中，只有能量项是正确的，用于统计的 sigma 和 E^2 并没有更新，所有需要使用其他工具，如 g_analysis 进行新的统计 （8）GMX.4 editconf editconf 是 GMX 最重要的程序之一。 它的主要功能是对系统结构进行编辑，同时它也把系统结构文件保存或者转换到不同的文件格式中，如 gro、g96、pdb 文件等。 在 分子动力学模拟中，通常给模拟系统添加一个模拟盒子（周期性的盒子，原子从这边出去了，就从那边进来，通俗吧）。 editconf 使用 “ -box ”、“ -d ” 和 “ -angle ” 等参数对模拟系统的盒子进行设定。在为系统设定盒子的时候，“ -box ” 和 “ -d ” 都把系统放置在盒子的中间，除非 editconf 明确使用另外一个参数 “ -noc ”。 (也就是 not center 的意思啦。) 使 用 “ -bt ” 参数，editconf 可以设定使用盒子的类型。editconf 支持以下几种盒子类型：triclinic （斜方体）、cubic（正方体）、dodecahedron (等边十二面体)， octahedron (等边八面体，两头被剪切的那种，即两个四面体放在一起，切去方向相反的尖尖的两头，同时保证所有的边长相等）。等边十二面体和等边八面体的体积分别是同“ 周期映像距离”的正方体体积的 0.71 和 0.77 倍，越接近球形，体积越小，计算代价越小。 当使用立方体、等边十二面体或者等边八面体时，“ -box ” 的参数可以只为一个实数值，该数值即为边长。“ -box ”参数也可以是三个实数，结果为以这三个实数为边长的长方体（如不指定 “ -angle ”，默认盒子矢量间夹角为90度）。 参数 “ -angle ” 只同参数 “ -box ” 和参数 “ -box ” 一起使用，以指定盒子矢量间的夹角；不能和参数 “ -d ”（用于指定系统原子到盒子边界的最小距离）一起使用。 如果使用 “ -n ” 或者 “ -ndef ” 参数为程序指定一个索引文件，则 editconf 可以选择系统中某一个组计算盒子大小或者几何中心，否则整个系统都会被考虑在内。 参数 “ -rotate” 用来旋转系统，按照给定数值在该坐标系上进行旋转。如 “ -rotate 0 30 0 ”，则表示将系统绕 Y 轴顺时间方向旋转 30 度。 参数 “ -princ ” 用来将系统（或者系统某一部分）固有坐标对齐到坐标轴上。这样做可以缩小盒子的体积，特别是当分子为一长条形状的时候。但是分子在模拟过程中会平移或者转动，所以使用的时候要特别小心。 扫 描参数 “ -scale ” 在其他参数被读入前可以扫描系统中各个原子的坐标，与其他参数配合使用可以对系统中原子坐标进行修改以得到不同的系统性质。如果 “ -density ” 一起使用，可以得到不同的系统密度（这样做会改变系统盒子大小，要特别小心！）。使用该参数时，如果输入为 .gro 文件，可能结果会不精确。 当 “ -scale ” 参数为单一 “ -1 ” 时，输出的系统结果为该方向上的镜面映像，当三个方向上的输入全部为 “ -1 ”（即 “-scale -1 -1 -1”）时，系统的结果机构为原来结构的坐标原点对称映像。 在程序输出时，可以只输出系统的某一个组，或者某一个部分。可以建立划分更为细致的索引文件，这样可以进行细致的选择。 系统结构的周期性可以在程序中进行粗略去除，但是必须保证结构文件中周期性盒子的信息绝对正确，因为 editconf 去除周期性的算法十分简单，只是将原子坐标直接减去盒子边长。 当 程序输出文件为 .pdb 文件时，可以使用 “ -bf ” 参数为输出文件添加 B-factor。B-factor 要在一定的格式的文件中读取，这个格式 要求 对教师党员的评价套管和固井爆破片与爆破装置仓库管理基本要求三甲医院都需要复审吗 文件第一行为文件中所含 B-factor 数值个数，文件第二行开始没行要有一个索引号，然后就是该索引后面的 B-factor 值。 B-factor 默认是一个残基一个数值，如果每一个原子有一直 B-factor 值，则要使用 “ -atom ” 参数。 如果使用 “ -legend ” 参数，那么程序在结构文件中会生成一列 CA 原子，这些原子带有 B-factor数字，并从最小 B-factor 值排列到最高，便于可视化软件显示。 如果使用 “ -mead ” 参数，那么程序将输出一个被静电统计软件 MEAD（玻松-波尔兹曼方程求解软件）使用的 pdb 文件（或者 pqr 文件）。使用这个参数有一个前提条件，就是输入文件必须是模拟打包文件（如 tpr 文件），因为这样的文件包含了力场参数。这个输出 pdb (pqr) 文件中，B-factor 数值列为该行原子的范德华半径，而 occupancy 列则为该原子的电量。“ -grasp ” 参数的作用很类似，但是 B-factor 列和 occupancy 列的数值互换。 另外一个十分的有用的参数就是 “ -label ” 参数，它可以为 pdb 文件加上一个链标记。如果一个文件里面不同残基属于不同肽链，那么这个参数可以帮这些残基指定肽链归属，这样不但可以帮助可视化，在建立模拟系统时也十分方便。 editconf 还可以修改结构文件的盒子类型，比如以下就是立方体盒子修改成一个等边八面体: +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ editconf -f -rotate 0 45 35.264 -bt octahedron -box -o +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 其中， 的数值为正方形变成乘于根号三除以二: a*sqrt(3)/2 。 程序输入文件： ------------- -f: 输入结构文件，文件格式可以是 gro、g96、pdb、tpr、tpb 或者 tpa。 -n: 索引文件。 -o: 输出文件，文件格式为 gro、g96 或者 pdb 。 -mead: 输出文件，MEAD程序的坐标文件，包含力场参数。 -bf: B-factor数值文件，文件格式见上文。 其他参数： --------- -w: 程序结束自动打开输出文件。 -ndef: 在索引文件的默认组中选择输出范围。 -bt: 指定盒子类型，支持的盒子类型有 triclinic（斜方体）、cubic（正方体）、dodecahedron（等边十二面体）或者 octahedron（等边八面体）。 -box: 指定盒子三个矢径，（a，b，c）。 -angle: 盒子矢径间夹角，默认都是 90 度（bc，ac，ab）。 -d: 溶质分子（蛋白，核酸等）到盒子边界的距离。 -c: 把溶质分子放在盒子中间，“-box” 和 “-d” 参数的默认参数。 -center: 为输出结构指定几何坐标中心，默认为 （0，0，0）。 -translate: 平移分子机构，比如 “-translate 0 0 2.5” 表示将分子结构向 Z 轴正方向平移 2.5 纳米。 -rotate: 旋转分子结构，见上文。 -princ: 对齐分子结构，见上文。 -scale: 扫描分子结构，“-scale 1 1 1” 表示在三个坐标方向上进行扫描。 -density: 指定输出结构文件的分子密度。 -pbc: 处理结构文件的周期性。 -grasp: 将原子的力场参数存储在输出文件中，见上文。 -rvdw: 给定默认原子的默认范德华半径。该参数可用于处理缺少力场参数的原子。 -sig56: 使用范德华能量最低点到原子中心距离的一半距离，而不使用简单使用范德华半径的一半。（特殊使用，多特殊？） -vdwread: 从默认 vdwradii.dat 文件中读取原子的范德华半径，而不使用力场提供的范德华半径。 -atom: 为各个原子添加 B-factor 数值。 -legend: 在结构文件中添加 B-factor 参考值，见上文。 -label: 为各个残基指定肽链归属。 注意：对于复杂的分子，最好使用 trjconv 去除结构周期性。 ：） (9)朗格文动力学 ( Langevin Dynamics ) 2009-06-19 14:10 朗格文动力学 ( Langevin Dynamics ) 是控制模拟系统能量的一种常用算法，在多种分子模拟软件中都可以看到。 分子模拟在一定的系综下进行，所以要保持系统状态不变，比如控制系统温度，压强等。由于计算机不是百分百精准，微小的误差在长时间的模拟过程中可能被不停积累和放大，于是需要不同的方法对系统进行不停调整。这些调整的方法有很多，比如 Berendersen, Nose-hoover, Langevin。 模拟系统中的原子受到起周围原子相互作用力，依据所受作用力根据牛顿第二定律运动。朗格文动力学的实现方法是给原子添加两个额外的作用力，即摩擦力和随机力。该摩擦力大小为原子速度乘于其质量再乘于一个摩擦因子（gamma)，其方向与原子速度相反。而随即力可以理解为来自溶液分子的随机相互作用等。这个两个力一起调节系统中各个原子的运动，以达到对整个系统能量的调控，即调控系统温度，压强等。 (10GMX.3 do_dssp 本程序调用第三发软件 DSSP，从轨迹文件中读出蛋白质的二级结构信息。要使用该程序，必先安装好 DSSP 软件（恩，用 google 或者 baidu 小搜一下）。安装 DSSP 完毕后(/home/user/soft/bin/dssp)，在 shell 中定义 DSSP 变量： setenv DSSP /home/user/soft/bin/dssp 如果使用 bash，则 export DSSP='/home/user/soft/bin/dssp'。可以该变量直接加到 shell 的配置文件中。 程序输出文件为 .xpm 文件，该文件其实是一个文本文件，可以一般文本编辑器打开。文件中用不同字符表示蛋白质每一个残基属于什么二级结构，并随时间变化，同时定义每个字符的颜色。.xpm 文件可视化有一点麻烦，Gromacs 的另外一个程序 xpm2ps 可以将其转化为 PostScript 文件，为 .eps 文件后缀，可以直接打开或直接使用到 Latex 文件中。 do_dssp可以统计属于某二级结构类型的残基数目，“-sss” 参数用于选择二级结构类型，“-sc“ 参数则把统计结果输出到 scount.xvg 文件中（其实这个文件包含了所有不同二级结构氨基酸残基数目，可以用 xmgrace 的 “-nxy” 参数打开）。本程序还可以计算每一个残基的溶液可及化面积(Solvent Accesible Surface, SAS)，以绝对平方埃（A^2)表示，或者最大 SAS 面积的分数表示（最大 SAS 定义为该一残基在甘氨酸链中的 SAS 面积）。Gromacs 另外一个程序 g_sas 也可以求取 SAS， 效率较高并且相对简单。 程序同时可以使用输出一个 ssdump.dat 文件。该文件包含了各个残基的二级结构信息，可以被另外一个程序 g_chi 使用，以分析残基直接二面角性质与二级结构直接的关系。（其实这个文件就是一个文本文件，里面用字符代表残基的二级结构，如 H 表示螺旋，B 表示折叠等。） 程序输入输出文件： --------------- -f: 轨迹文件，xtc、trr等。 -s : 模拟系统 tpr 文件。 -n: 索引文件。 -ssdump: 输出 ssdump.dat 文件。 -map: 程序输出 xmp 文件的原色库文件，如无则默认输出。（让其默认挺好的啦） -o: 输出文件，xmp 文件。各个残基属于某二级结构信息并随模拟时间变化。 -sc: 统计某二级结构残基数目，输出文件 scount.xvg。 -a: 各个残基的溶液可及化面积，xmp文件。 -ta: 总可及化面积，包括疏水和亲水面积，xvg 文件。 -aa： 平均可及化面积，xvg 文件。 其他参数： -------- -b: 指定可是统计时间帧。 -e: 指定结束统计时间帧。（“-b” 和 “-e” 可以用来读取轨迹中某一时间断的信息。）-dt: 指定读取帧的时间间隔。 -tu: 设定输出时间单位。 -w: 统计结束，自动打开输出文件。 -xvgr: 在输出 xvg 文件中添加其他信息，坐标标题等。 -sss: 设定统计二级结构类型。 (11GMX.2 anadock anadock是一个分析分子对接（docking）软件 Autodock 计算结果的程序，其作用是根据距离或者 RMSD 对分子对接结果的分子结构进行归簇，然后依据分子对接能量和自由能，并进行统计和输出。（本人没有用过 AutoDock，只用过 Dock，故不作任何评论。） 对分子结构进行归簇，Gromacs 也提供 g_cluster 程序，故可以使用该程序对结构进行提前归簇在计算能量，同样可以达到目的。 程序输入输出文件： ----------------------- -f: 输入文件，PDB 文件格式。 -ox: 归簇输出文件，PDB文件格式。 -od: 能量输出文件，为 xmgrace 的 xvg 文件格式。 -of: 自由能输出文件，xvg 文件格式。 -g: 程序记录文件。 其他参数： ------------ -xvgr: 在 xvg 输出文件中添加文字 说明 关于失联党员情况说明岗位说明总经理岗位说明书会计岗位说明书行政主管岗位说明书 ，坐标名称等。 -free: 使用 AutoDock 估算的自由能对结构进行结构排序。 -rms: 使用 RMS 进行归簇。 -cutoff: 归簇尺度，偏离大于给出数值时认为不同簇，默认为 0.2 nm。 （13Restraint & Constraint 在分子模拟中，这两个词都表示限制，但是意义不同。通俗点讲，受到 Restraint 的分子或原子原则上可以运动；受到 Constraint 的分子原子原则上不能运动。 Restraint 通常的实现方法是给目标施加一个弹性势能，当目标偏离平衡位置时，弹性势能使目标恢复到平衡位置。Constraint 则直接使目标固定住，则施加无限大的弹性势能。 在 MD 软件如Gromacs中，Restraint 定义的弹性势能直接施加到指定原子上，并在 pores.itp 文件中列出。当 MDP （Molecular Dynamics Parameters) 文件使用关键字 “define = -DPOSRES” 时，引用该 pores.itp 文件。 Constraint 的实现则在MPD文件中定义 ”freezegrps “ 关键字，把需要进行 Constraint 的分子或原子直接添加到关键字中。 BTW: GMX 中 Restraint 和 Constraint 都具有方向性，即限制可以在某一个方向上实现。其他模拟软件雷同。 （14GMX.1 参数概述 根据GMX手册预告，所有GMX的程序都有6个标准参数，其中有的是隐藏的（参数太多，吓死人）。此六个参数及其作用罗列如下： -h: 在屏幕上打印该GMX程序的帮助信息，并退出程序。就是不使用这个参数，任何GMX程序运行时，都会打印很多信息，包括GMX的软件消息，简短的 协议 离婚协议模板下载合伙人协议 下载渠道分销协议免费下载敬业协议下载授课协议下载 ，还有各个参数的输入等。 -hidden: 打印隐藏的命令参数。 -quiet: 打印尽可能少的信息，所有参数输出被忽略，所以很quiet。 -man: 以指定文件格式，打印程序的手册到指定文件中。支持的文件格式有html, tex, nroff, acsii, completion, py 和 xml. （真怀疑开发人员怎么不被这么多文件格式给烦死，牛人。） -debug: 输出程序的调试信息。这个对终端用户不是很有用，输出文件一般为core加一串数字。 -nice: 设置程序的优先级。一般忽略，相信聪明可爱的计算机。 以下为GMX一些其他介绍： 1、GMX编译过程中，如果使用Motif或者Lesstif图形包，那么可以使用 -X 参数得到一个编辑各个参数输入的图形界面。如果没有使用以上图形包，使用该参数会报错。 2、在SGI-IRIX系统上编译的GMX，可以使用 -npri 参数设定无阻塞（non-blocking）优先级。 3、可选参数的如无指定文件，那么就不会有输出文件（不然垃圾文件太多）。但是必选参数的输出文件一定会输出，如果指定按照默认文件明输出。 4、所有的GMX程序都接受无文件后缀名的文件，并自行补全文件名。如果不指明输入文件命，程序会自动在当前目录下搜索默认文件名。如果有多个同名不同文件类型的输入，那么第一个支持类型文件会被使用。（相信很多人一般不至于懒惰到这个程度，连个文件名都不给人家程序。） 5、几乎所有GMX程序（不包括 mdrun、nmrun 和 eneconv）都会先检查输入参数是否合法（重复参数等），如出错程序自动退出。 6、枚举参数（enum）必须为指定类型之一。 7、矢量输入必须文一个或者三个数值，如果指定一个数值，那么其他两个与指定数值一致。 8、多数GMX程序都可以使用 -tu 参数来指定时间单位，如 ps、fs、ns等等。这个多数应用与分析输出，成图等。 9、所有GMX程序可以直接读取 g-zipped 压缩文件，xtc、trr 和 trj 等压缩轨迹文件可能会报错，建议只打包不压缩。 10、大多数GMX程序可以处理比输入文件中原子数少的轨迹文件（轨迹文件中原子数目少于输入文件），但是必须保重轨迹文件中的原子处于输入文件的开头（即轨迹文件中没有的原子在输入文件尾端）。（这个貌似很容易产生歧义，语文没学好，见谅！） 造成这个的原因是GMX系统中，输入文件中所有原子按照序号进行，只有一个数字，没有名称。 （15进行模拟的基本步骤行 以下是进行分子动力学模拟的最基本步骤，来自Gromacs维基。具体的步骤可能会因为模拟目标的不同而存在差异，这里只做基本的概述。 十分清楚要什么样的模拟系统性质和现象，明白模拟的目标。 选择适当的工具以实现模拟目标。这需要熟知其他研究人员对相似模拟系统的模拟方法等资料，要好好读别人的论文，主要包括： -模拟使用的软件，当然还要考虑使用的力场，有时因为力场关系，需要更改使用的软件。 -力场描述了系统中原子的性质和相互作用，选择一个适合自己研究目标的力场至关重要，还是需要多看论文。根据sensenbobo血泪史和别人的经验，模拟核酸的话，可以采用amber力场，水分子可以使用TIP3P；模拟蛋白质的话，可以采用OPLS-AA，水分子使用TIP4P（重申一下，我本人没有跑过测试）。 获取模拟分子的坐标文件（PDB），或者自己生成一个。自己生成很好啊，强烈鼓励，但是要十分明白有一定的危险性，要明白每一步模拟系统的现象是为什么。 根据模拟软件，生成模拟系统初步坐标文件。 获取或者生成系统的拓扑文件，比如Gromacs的pdb2gmx命令，或者使用网络服务器PRODRG（google一下这几个字母），amber的leap命令和antechamber，NAMD的psfgen或者VMD。请参考相关手册。 给模拟系统添加一个盒子（盒子的翻译方法我十分不喜欢，箱子显得太大，边界显得太难理解，但是盒子听起来想到蛋糕，我们是严肃的科学家。娃哈哈。）在盒子中添加水分子，添加离子如NA和CL之类（添油加醋，生活就是这样进行的）。添加任何东西到你的系统中，就要相应的修改一下拓扑文件，以保持与系统对应，不然会出错哦。 进行能量最优话。这样做是消去系统中原子之间距离太近造成的高能量。必要的时候，分子可以现在真空中能量最优话，然后再做溶剂环境中能量最优话。 选择恰当的模拟参数模拟平衡模拟系统。需要慎重考虑力场的作用，必要的时候保持系统体积温度不变（NVT系综）进行位置限制性平衡系统，然后在保持系统温度压强不变（NPT系综）修正系统的密度。然后再选择正确的系综（NPT，NVT，NET，BT？？）进行数据收集。（真烦！） 跑足够长时间的模拟，使系统达到理想的平衡。使用模拟软件的工具或者其他软件观察系统的性质，如温度，能量瞪瞪。 使用适当的参数进行数据采集模拟，注意模拟衔接过程中各个系统性质保持不变，如果原子速度之类。这里还要考虑力场的作用，同时要考虑怎么测量描述模拟目标。 模拟足够长的时间使模拟目标明显出现（可能也出现跟你期望相反的结果，出现率应该百分之五十，加上个人感情因素，可能会达到百分之八十。保重！）。 分析模拟结果，解释得到结果。要好好利用视图工具，同时要十分注意漂亮的三维图说明不了什么问题。图表，比较，统计，一个都不能少。 就这样吧。 （16Gromacs使用－新手入门篇】 其实gromacs的思路比较简单，大体就是以下的流程,在此先对em进行一下 总结 初级经济法重点总结下载党员个人总结TXt高中句型全总结.doc高中句型全总结.doc理论力学知识点总结pdf ： 1. convert the pdb-file to a gromacs structure file(.gro / .pdb) and a gromacs topology file(.top) 这一步需要gromacs的预处理程序pdb2gmx 标准的输入命令应该是： ＃ pdb2gmx -ignh -f speptide.pdb -p speptide.top -o speptide.gro 这个命令之后有一个force field的选择 系统给出了11种力场，具体以后实际计算中如何选择力场还是应该多参阅文献和gromacs说明书的。 这步之后可以在屏幕看到： --------- PLEASE NOTE ------------ You have succesfully generated a topology from: speptide.pdb. The G43a1 force field and the spc water model are used. Note that the default mechanism for selecting a force fields has changed, starting from GROMACS version 3.2.0 --------- ETON ESAELP ------------至此，完成了从pdb到gro和top的转变。另外，我也看过有的作这一步命令时候并不进行top文件的转换，因为在genbox（对盒子加水溶剂化时还要对top文件进行修改） 2. solvate the peptide in water 这一步需要用到editconf和genbox两个程序，命令的标准书写应该是： # editconf -bt cubic -f speptide.gro -o speptide_box.gro -c -d 0.5 可以在屏幕上看到 Read 191 atoms Volume: 8.17315 nm^3, corresponds to roughly 3600 electrons No velocities found system size : 1.774 3.372 1.367 (nm) diameter : 3.517 (nm) center : 2.650 1.453 2.417 (nm) box vectors : 1.774 3.372 1.366 (nm) box angles : 90.00 90.00 90.00 (degrees) box volume : 8.17 (nm^3) shift : -0.391 0.806 -0.158 (nm) new center : 2.259 2.259 2.259 (nm) new box vectors : 4.517 4.517 4.517 (nm) new box angles : 90.00 90.00 90.00 (degrees) new box volume : 92.18 (nm^3) 因为程序的默认，所以这一步命令可以写作： # editconf –f speptide –o –d 0.5 这里－f之后的sp文件是默认gro文件或者pdb文件，－o后面生成的output file是放入盒子的体系的gro文件 按照这样的格式，输出的文件被默认命名为－out.gro 体系放入盒子之后在盒子中注入水 # genbox -cp speptide_box.gro -cs -p speptide.top -o speptide_water.gro 也可简写为 ： # genbox –cp out –cs –p sprptide –o b4em Screen-- Output configuration contains 9035 atoms in 2967 residues Volume : 92.1813 (nm^3) Density : 994.449 (g/l) Number of SOL molecules: 2948 之后看加水之后的top文件 用命令 # tail speptide.top 可以看到 [ system ] ; Name Protein in water [ molecules ] ; Compound #mols Protein 1 SOL 2948 其中SOL后面的数字2948就是genbox加入盒子的水分子的数目。 Editconf程序的另一个用途是讲gro文件转化回pdb这时可以讲speptide_water.gro转化回pdb观察 ＃ editconf –f speptide_water.gro –o speptide.pdb 拖回本机 用spbdv或者vmd观察 。 The next step is to generate index file . in the tutor of gromacs , we are told that there are a set of index groups to select, unfortunately, I didi not find them ,so I have to use make_ndx to generate one . # make_ndx –f b4em 3. perform an enery minimization of the peptide in solvent Now the simulation system is almost ready. Before we can start the dynamics, we must perform an energy minimization, to alleviate any bad contacts (atoms overlapping such that a significant repulsion would result, causing numerical problems in the simulation) that might be present in the system. # grompp –v –f em.mdp –c speptide_water.gro –p speptide.top –o speptide_em.tpr Or use this command for short: # grompp –v –f em –c b4em –p speptide –o em After this command ,the bad contacs have been removed .so we can do enery minimization now. # mdrun –v –s speptide_em.tpr –o speptide_em.trr –c after_em.gro –g emlog.log Or # mdrun –v –s em –o em –c after_em –g emlog From screen ,we can see these ------ Steepest Descents: Tolerance (Fmax) = 2.00000e+03 Number of steps = 100 Step= 0, Dmax= 1.0e-02 nm, Epot= -8.02562e+04 Fmax= 2.70468e+04, atom= 6792 Step= 1, Dmax= 1.0e-02 nm, Epot= -8.69641e+04 Fmax= 1.24531e+04, atom= 6657 Step= 2, Dmax= 1.2e-02 nm, Epot= -9.48067e+04 Fmax= 5.54072e+03, atom= 9009 Step= 3, Dmax= 1.4e-02 nm, Epot= -1.02743e+05 Fmax= 2.86957e+03, atom= 8766 Step= 4, Dmax= 1.7e-02 nm, Epot= -1.08941e+05 Fmax= 1.53075e+04, atom= 150 Step= 5, Dmax= 2.1e-02 nm, Epot= -1.10083e+05 Fmax= 1.31653e+04, atom= 150 Step= 6, Dmax= 2.5e-02 nm, Epot= -1.10564e+05 Fmax= 2.16454e+04, atom= 150 Step= 7, Dmax= 3.0e-02 nm, Epot= -1.11528e+05 Fmax= 1.82398e+04, atom= 150 Step= 9, Dmax= 1.8e-02 nm, Epot= -1.12767e+05 Fmax= 3.87548e+03, atom= 120 Step= 11, Dmax= 1.1e-02 nm, Epot= -1.13598e+05 Fmax= 1.23013e+04, atom= 120 Step= 12, Dmax= 1.3e-02 nm, Epot= -1.14464e+05 Fmax= 5.37656e+03, atom= 120 Step= 13, Dmax= 1.5e-02 nm, Epot= -1.14787e+05 Fmax= 1.74494e+04, atom= 120 Step= 14, Dmax= 1.9e-02 nm, Epot= -1.15874e+05 Fmax= 8.35820e+03, atom= 120 Step= 16, Dmax= 1.1e-02 nm, Epot= -1.16420e+05 Fmax= 6.65421e+03, atom= 120 Step= 17, Dmax= 1.3e-02 nm, Epot= -1.16807e+05 Fmax= 1.08106e+04, atom= 120 Step= 18, Dmax= 1.6e-02 nm, Epot= -1.17248e+05 Fmax= 1.12834e+04, atom= 120 Step= 19, Dmax= 1.9e-02 nm, Epot= -1.17425e+05 Fmax= 1.45145e+04, atom= 120 Step= 20, Dmax= 2.3e-02 nm, Epot= -1.17590e+05 Fmax= 1.67089e+04, atom= 120 Step= 22, Dmax= 1.4e-02 nm, Epot= -1.18845e+05 Fmax= 2.65542e+03, atom= 121 Step= 23, Dmax= 1.7e-02 nm, Epot= -1.19198e+05 Fmax= 2.39072e+04, atom= 120 Step= 24, Dmax= 2.0e-02 nm, Epot= -1.20588e+05 Fmax= 6.17310e+03, atom= 120 Step= 26, Dmax= 1.2e-02 nm, Epot= -1.20901e+05 Fmax= 1.04224e+04, atom= 120 Step= 27, Dmax= 1.4e-02 nm, Epot= -1.21229e+05 Fmax= 8.89909e+03, atom= 120 Step= 28, Dmax= 1.7e-02 nm, Epot= -1.21346e+05 Fmax= 1.51125e+04, atom