请选择 进入手机版 | 继续访问电脑版

计算化学公社

 找回密码
 现在注册!
查看: 13958|回复: 15

[VMD] 在VMD中计算RMSD衡量两个结构间的差异以及叠合两个结构

[复制链接]

2万

帖子

25

威望

3万

eV
积分
55897

管理员

公社社长

发表于 2015-5-11 01:04:06 | 显示全部楼层 |阅读模式
在VMD中计算RMSD衡量两个结构间的差异以及叠合两个结构

文/Sobereva @北京科音
First release: 2015-May-11  Last update: 2020-Feb-25


1 原理

RMSD(Root mean square displacement/deviation)搞分子模拟的人都很熟悉,但是搞量子化学的人可能好多都不知道。有人问怎么衡量激发态结构相对于基态的改变,虽然可以挨个比较键长键角的变化,但是若想用一个值衡量整体的结构改变程度,计算RMSD是最恰当也是最省事的,这里就简单说一下怎么算。对任意两个结构间都可以计算RMSD,比如可以衡量两种极小点构型的差异、电子激发/电离/外加电场前后几何结构的整体改变程度、形成复合物后单体结构的变化等等。

RMSD定义为
1.png

其中i循环所有原子,x_i和x_i'分别是第i个原子在第一个结构和第二个结构中的x坐标,y、z类似。

计算两个结构间的RMSD之前必须先进行叠合(align),相当于令一个结构平移和旋转来最小化两个结构之间的RMSD,不做这一步的话计算出的RMSD是不能衡量两个结构内部的差异的。

能计算RMSD、做叠合的程序不计其数,下面用免费的VMD程序(http://www.ks.uiuc.edu/Research/vmd/)的1.9.3版来演示,而且还将体现如何通过叠加图直观地展现两个结构之间的差异。


2 实例:中性和阳离子状态下的丙烷

本例演示如何计算丙烷在中性状态下和在+1阳离子状态下优化的结构之间的RMSD值,并且绘制出叠合后的图。

首先把中性的丙烷用量子化学程序优化,把优化后的结构保存成VMD能认的一种格式,比如pdb、mol2、xyz等。这里我们用很常用的xyz格式,见《谈谈记录化学体系结构的xyz文件》(http://sobereva.com/477)。我们创建一个名为neutral.xyz的文本文件,把优化后的中性丙烷坐标拷进去,并在第一行写上原子数,第二行是注释行(内容随便写)。具体内容如下
    11
Neutral
6        0.000000    1.277193   -0.259856
1        0.884678    1.321887   -0.907445
1       -0.884678    1.321887   -0.907445
1        0.000000    2.176504    0.367166
6        0.000000    0.000000    0.586502
1        0.877607    0.000000    1.247354
1       -0.877607    0.000000    1.247354
6        0.000000   -1.277193   -0.259856
1        0.000000   -2.176504    0.367166
1       -0.884678   -1.321887   -0.907445
1        0.884678   -1.321887   -0.907445


类似地也优化丙烷阳离子,把结构存为cation.xyz,内容如下
    11
Ionized
6       1.249671   -0.343789    0.000081
1       1.315464   -0.947246    0.909976
1       1.314553   -0.947793   -0.909563
1       2.180924    0.290109   -0.000844
6       0.192643    0.673245    0.000031
1       0.022018    1.237089    0.919487
1       0.022792    1.237187   -0.919539
6      -1.442724   -0.286303   -0.000053
1      -2.107101    0.577835   -0.000705
1      -1.372956   -0.853783   -0.925484
1      -1.373233   -0.852316    0.926316


PS:如果你用的是Gaussian,把优化任务最后一帧的结构转化为xyz格式最简单的方法是用Multiwfn(http://sobereva.com/multiwfn)。把Multiwfn的settings.ini文件里的iloadGaugeom设为1,然后用Multiwfn载入Gaussian优化任务的输出文件,选择主功能100的子功能2,就可以看见导出xyz文件的选项,选择即可。

启动VMD,把neutral.xyz和cation.xyz都拖到VMD Main窗口里载入。然后选择Extensions - Analysis - RMSD Calculator,把左上角的文本框的内容改为all,取消选择Backbond only复选框,然后点击Align,在图形窗口你会发现两个结构已经叠合了,重叠度较高。然后再点RMSD按钮,在RMSD Tool窗口下方显示的Total RMSD就是两个结构间的RMSD值。当前窗口看到的情况如下所示,可见RMSD是0.093埃。
2.png


做叠合和计算RMSD的时候,我们可以只对指定部分进行。比如我们不输入all,而是输入比如serial 2 to 4 8,那么点Align按钮时只会根据2、3、4、8号原子来做align,点击RMSD按钮时也只会计算这些原子间的RMSD。如果你对VMD的选择语句不熟悉,建议参看《VMD里原子选择语句的语法和例子》(http://sobereva.com/504)。

VMD还可以把两帧结构用不同颜色绘制出来便于直观比较。做法是选Graphics-Representation,在Selected Molecule里先选择对应中性的体系,把Coloring Method设为Color ID,在旁边选择13 Mauve,并把Thickness设为4使细线显示得更粗。然后在Selected Molecule里选择对应阴离子的体系,以类似操作把颜色设为蓝色并加粗细线。最后在文本窗口输入color Display Background white命令把背景改为白色,此时在图形窗口会看到下图,可见结构差异体现得十分清晰直观

3.png



3 实例:不同的取代苯体系

有时候我们想对两个主体结构相同,但某些部分不同的两个分子的局部区域计算RMSD来衡量局部结构的差异。本例就用苯酚以及间硝基苯腈为例进行说明,二者在B3LYP/6-31G*优化后的结构如下

4.png


二者的差异在于苯环上连的基团。本例我们要让两个分子的苯环部分相互叠合(只需考虑碳原子)并计算RMSD。为了能够正确叠合,一定要注意被叠合的部分要满足两个条件:(1)原子数相同 (2)原子顺序一致。由上图可见,这两个体系中的苯环上的碳都是按连接顺序排列的,所以可以直接做叠合和计算RMSD。如果其中一个体系的苯环上的碳是比如1,4,2,6,5,3这样交错排列的就没法直接做叠合了,而需要在叠合之前自行编辑结构文件,调整原子顺序使得条件(2)被满足。

还是先载入这两个结构到VMD里并进入RMSD Calculator,左上角的文本框输入serial 1 to 6(对应两个分子的苯环上的碳原子序号),取消Backbone only复选框,然后点Align。按理说此时两个分子的苯环应该被叠合了、之后再点击RMSD就完事了,然而当前的情况比较特殊,图形窗口里目前看到如下情况

5.png


显然还没有真正叠合好。这是VMD用的算法的缺陷,类似于优化落到了局部极小点而非全局最小点的意味。解决办法就是手动旋转其中一个分子,令其苯环部分的各个原子与另一个分子的苯环部分相应序号的原子整体比较接近,类似于提供一个更好的初猜。调节结构具体的做法是选Mouse - Move - Molecule,之后按住alt键并用鼠标左键拖动某原子就可以平移整个体系,按住alt+shift并用鼠标左键拖动某原子可以旋转整个分子,按住alt+shift并用鼠标右键拖动某原子可以令整个分子沿屏幕旋转。等调节得差不多了,再次按Align按钮,就会发现两个分子的苯环部分确实很好地叠合在了一起,然后再点击RMSD按钮,此时情况如下所示。

6.png


可见两个分子的苯环碳原子间的RMSD仅为0.013埃,体现出苯环非常刚,取代基对这部分结构的影响可以忽略不计。

评分

参与人数 2eV +10 收起 理由
yesheng + 5 好物!
kunkun + 5 とてもいい!

查看全部评分

北京科音自然科学研究中心http://www.keinsci.com  致力于计算化学的发展和传播,长期开办最高水准的各种量子化学、分子动力学、波函数分析与Multiwfn程序等主题的培训,是提升计算化学研究水平的最佳选择。欢迎加入“北京科音”微信公众号获取培训最新消息和计算化学资讯!培训相关信息见《北京科音办的培训班FAQ》(http://bbs.keinsci.com/thread-5098-1-1.html)。
欢迎加入人气最高、水准最高的综合性理论与计算化学交流QQ群:思想家公社QQ群1号:18616395,2号:466017436。合计约6000人。两个群讨论范畴相同,可加入任意其一但不可都加入,申请信息必须注明具体研究方向,否则一定会被拒绝加入。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(最强大的量子化学波函数分析程序)
ResearchGate:https://www.researchgate.net/profile/Tian_Lu

Money and papers are rubbish, get a real life!

186

帖子

0

威望

2090

eV
积分
2276

Level 5 (御坂)

发表于 2016-3-22 10:56:56 | 显示全部楼层
每次运算完一次RMSD计算后,都需要重启一次VMD,不然报错,怎样清楚上一次作业的痕迹?

129

帖子

0

威望

946

eV
积分
1075

Level 4 (黑子)

发表于 2017-6-20 16:23:36 | 显示全部楼层
老师我有两个问题:
1.在这里讲的是两帧构型的的计算,那如果我想做更多个,比如我在.xyz文件中保存100个分子,根据上述描述进行同样动作的操作,VMD可不可以实现对于多个构型RMSD值的计算?
2.老师您提到,在计算RMSD值之前要先对分子进行align,也就是通过平移翻转尽可能和第一帧构型接近,那对于这一部分程序内部是通过什么算法实现的,有没有计算RMSD的align的代码,我想学习一下?
谢谢老师~

105

帖子

1

威望

1115

eV
积分
1240

Level 4 (黑子)

发表于 2017-6-20 22:40:07 | 显示全部楼层
jackie 发表于 2017-6-20 16:23
老师我有两个问题:
1.在这里讲的是两帧构型的的计算,那如果我想做更多个,比如我在.xyz文件中保存100个 ...

关于第二个问题 可以参考https://github.com/charnley/rmsd 里面也有算法的介绍
第一个问题 如果你会python的话 调用上面的那个程序也很方便

129

帖子

0

威望

946

eV
积分
1075

Level 4 (黑子)

发表于 2017-6-21 10:26:32 | 显示全部楼层
KiritsuguPapa 发表于 2017-6-20 22:40
关于第二个问题 可以参考https://github.com/charnley/rmsd 里面也有算法的介绍
第一个问题 如果你会pyt ...

So cool!!!非常感谢

2万

帖子

25

威望

3万

eV
积分
55897

管理员

公社社长

 楼主| 发表于 2017-6-21 13:48:26 | 显示全部楼层
jackie 发表于 2017-6-20 16:23
老师我有两个问题:
1.在这里讲的是两帧构型的的计算,那如果我想做更多个,比如我在.xyz文件中保存100个 ...


RMSD Trajectory Tool里,点了align,所有帧(xyz里每个构型)都会相对于第一帧(默认情况下如此)做align,然后可以直接把所有帧的RMSD都算出来并导出
北京科音自然科学研究中心http://www.keinsci.com  致力于计算化学的发展和传播,长期开办最高水准的各种量子化学、分子动力学、波函数分析与Multiwfn程序等主题的培训,是提升计算化学研究水平的最佳选择。欢迎加入“北京科音”微信公众号获取培训最新消息和计算化学资讯!培训相关信息见《北京科音办的培训班FAQ》(http://bbs.keinsci.com/thread-5098-1-1.html)。
欢迎加入人气最高、水准最高的综合性理论与计算化学交流QQ群:思想家公社QQ群1号:18616395,2号:466017436。合计约6000人。两个群讨论范畴相同,可加入任意其一但不可都加入,申请信息必须注明具体研究方向,否则一定会被拒绝加入。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(最强大的量子化学波函数分析程序)
ResearchGate:https://www.researchgate.net/profile/Tian_Lu

Money and papers are rubbish, get a real life!

129

帖子

0

威望

946

eV
积分
1075

Level 4 (黑子)

发表于 2017-6-23 10:41:03 | 显示全部楼层
sobereva 发表于 2017-6-21 13:48
RMSD Trajectory Tool里,点了align,所有帧(xyz里每个构型)都会相对于第一帧(默认情况下如此)做al ...


老师 为什么我把好多帧的蛋白构型文件用VMD打开通过Graphics-Representation-Drawing method选项没法变成cartoon构型?文件已上传

2万

帖子

25

威望

3万

eV
积分
55897

管理员

公社社长

 楼主| 发表于 2017-6-23 15:50:30 | 显示全部楼层
jackie 发表于 2017-6-23 10:41
老师 为什么我把好多帧的蛋白构型文件用VMD打开通过Graphics-Representation-Drawing method选项没法变 ...

上传大文件前一定要压缩以节约论坛空间、加快下载速度
xyz文件里没有残基、原子名信息,不足以令VMD通过stride程序判断二级结构
北京科音自然科学研究中心http://www.keinsci.com  致力于计算化学的发展和传播,长期开办最高水准的各种量子化学、分子动力学、波函数分析与Multiwfn程序等主题的培训,是提升计算化学研究水平的最佳选择。欢迎加入“北京科音”微信公众号获取培训最新消息和计算化学资讯!培训相关信息见《北京科音办的培训班FAQ》(http://bbs.keinsci.com/thread-5098-1-1.html)。
欢迎加入人气最高、水准最高的综合性理论与计算化学交流QQ群:思想家公社QQ群1号:18616395,2号:466017436。合计约6000人。两个群讨论范畴相同,可加入任意其一但不可都加入,申请信息必须注明具体研究方向,否则一定会被拒绝加入。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(最强大的量子化学波函数分析程序)
ResearchGate:https://www.researchgate.net/profile/Tian_Lu

Money and papers are rubbish, get a real life!

50

帖子

0

威望

550

eV
积分
600

Level 4 (黑子)

发表于 2018-8-6 16:23:22 | 显示全部楼层
sobereva老师,打开VMD后,Extensions-Analysis目标下面有RMSD calculator,RMSD Trajectory Tool,RMSD visualizer Tool三个选项,这三个选项计算RMSD有什么区别吗?

2万

帖子

25

威望

3万

eV
积分
55897

管理员

公社社长

 楼主| 发表于 2018-8-7 04:12:35 | 显示全部楼层
陈丢丢 发表于 2018-8-6 16:23
sobereva老师,打开VMD后,Extensions-Analysis目标下面有RMSD calculator,RMSD Trajectory Tool,RMSD vi ...


原理都一样,只不过具体用处不同。一般只用得到RMSD Trajectory Tool
北京科音自然科学研究中心http://www.keinsci.com  致力于计算化学的发展和传播,长期开办最高水准的各种量子化学、分子动力学、波函数分析与Multiwfn程序等主题的培训,是提升计算化学研究水平的最佳选择。欢迎加入“北京科音”微信公众号获取培训最新消息和计算化学资讯!培训相关信息见《北京科音办的培训班FAQ》(http://bbs.keinsci.com/thread-5098-1-1.html)。
欢迎加入人气最高、水准最高的综合性理论与计算化学交流QQ群:思想家公社QQ群1号:18616395,2号:466017436。合计约6000人。两个群讨论范畴相同,可加入任意其一但不可都加入,申请信息必须注明具体研究方向,否则一定会被拒绝加入。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(最强大的量子化学波函数分析程序)
ResearchGate:https://www.researchgate.net/profile/Tian_Lu

Money and papers are rubbish, get a real life!

50

帖子

0

威望

550

eV
积分
600

Level 4 (黑子)

发表于 2018-8-7 09:53:37 | 显示全部楼层
sobereva 发表于 2018-8-7 04:12
原理都一样,只不过具体用处不同。一般只用得到RMSD Trajectory Tool

好的,谢谢sobereva老师!

2万

帖子

25

威望

3万

eV
积分
55897

管理员

公社社长

 楼主| 发表于 2020-2-25 07:54:59 | 显示全部楼层
大幅重写了本文,使用法更为灵活,并且加入了对原子数不同的两个体系做叠合并计算RMSD的说明
北京科音自然科学研究中心http://www.keinsci.com  致力于计算化学的发展和传播,长期开办最高水准的各种量子化学、分子动力学、波函数分析与Multiwfn程序等主题的培训,是提升计算化学研究水平的最佳选择。欢迎加入“北京科音”微信公众号获取培训最新消息和计算化学资讯!培训相关信息见《北京科音办的培训班FAQ》(http://bbs.keinsci.com/thread-5098-1-1.html)。
欢迎加入人气最高、水准最高的综合性理论与计算化学交流QQ群:思想家公社QQ群1号:18616395,2号:466017436。合计约6000人。两个群讨论范畴相同,可加入任意其一但不可都加入,申请信息必须注明具体研究方向,否则一定会被拒绝加入。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(最强大的量子化学波函数分析程序)
ResearchGate:https://www.researchgate.net/profile/Tian_Lu

Money and papers are rubbish, get a real life!

26

帖子

0

威望

137

eV
积分
163

Level 3 能力者

发表于 2020-5-19 11:55:40 | 显示全部楼层
你好老师,我想问一下,我有一个目标单体和一个已知三聚体结构,想通过叠加做出目标的三聚体结构,怎么操作呢。

2万

帖子

25

威望

3万

eV
积分
55897

管理员

公社社长

 楼主| 发表于 2020-5-21 08:01:38 | 显示全部楼层
小土豆 发表于 2020-5-19 11:55
你好老师,我想问一下,我有一个目标单体和一个已知三聚体结构,想通过叠加做出目标的三聚体结构,怎么操作 ...

我不知道你说的“已知三聚体结构”具体指什么,里面的分子和当前单体有什么差别
如果里面的分子没法和当前单体相叠合,不可能靠叠加得到当前单体的三聚体结构
北京科音自然科学研究中心http://www.keinsci.com  致力于计算化学的发展和传播,长期开办最高水准的各种量子化学、分子动力学、波函数分析与Multiwfn程序等主题的培训,是提升计算化学研究水平的最佳选择。欢迎加入“北京科音”微信公众号获取培训最新消息和计算化学资讯!培训相关信息见《北京科音办的培训班FAQ》(http://bbs.keinsci.com/thread-5098-1-1.html)。
欢迎加入人气最高、水准最高的综合性理论与计算化学交流QQ群:思想家公社QQ群1号:18616395,2号:466017436。合计约6000人。两个群讨论范畴相同,可加入任意其一但不可都加入,申请信息必须注明具体研究方向,否则一定会被拒绝加入。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(最强大的量子化学波函数分析程序)
ResearchGate:https://www.researchgate.net/profile/Tian_Lu

Money and papers are rubbish, get a real life!

149

帖子

0

威望

682

eV
积分
831

Level 4 (黑子)

发表于 2020-6-23 16:15:01 | 显示全部楼层
老师,我想引用一下RMSD计算公式,但也没找到您写的这个,有参考可以提供一下吗?
您需要登录后才可以回帖 登录 | 现在注册!

本版积分规则

手机版|北京科音自然科学研究中心|京公网安备 11010502035419号|计算化学公社 — 北京科音旗下高水平计算化学交流论坛 ( 京ICP备14038949-1号 )

GMT+8, 2020-7-4 23:52 , Processed in 0.180899 second(s), 27 queries .

快速回复 返回顶部 返回列表