计算化学公社

 找回密码 Forget password
 注册 Register
Views: 53641|回复 Reply: 35
打印 Print 上一主题 Last thread 下一主题 Next thread

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

  [复制链接 Copy URL]

6万

帖子

99

威望

5万

eV
积分
120134

管理员

公社社长

在VMD中计算RMSD衡量两个结构间的差异以及叠合两个结构
Calculating RMSD in VMD to measure the difference between two structures and to superimpose two structures

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


1 原理

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

RMSD定义为

其中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埃。


做叠合和计算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 实例:不同的取代苯体系

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



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

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



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



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

评分 Rate

参与人数
Participants 4
eV +20 收起 理由
Reason
人间划水艺术家 + 5 好萌好萌好萌!
邓苏微 + 5 赞!
yesheng + 5 好物!
kunkun + 5 とてもいい!

查看全部评分 View all ratings

北京科音自然科学研究中心http://www.keinsci.com)致力于计算化学的发展和传播,长期开办极高质量的各种计算化学类培训:初级量子化学培训班中级量子化学培训班高级量子化学培训班量子化学波函数分析与Multiwfn程序培训班分子动力学与GROMACS培训班CP2K第一性原理计算培训班,内容介绍以及往届资料购买请点击相应链接查看。这些培训是计算化学从零快速入门以及进一步全面系统性提升研究水平的高速路!培训各种常见问题见《北京科音办的培训班FAQ》
欢迎加入北京科音微信公众号获取北京科音培训的最新消息,并避免错过网上有价值的计算化学文章!
欢迎加入人气极高、专业性特别强的理论与计算化学综合交流群思想家公社QQ群(群号见此链接),合计达一万多人。北京科音培训班的学员在群中可申请VIP头衔,提问将得到群主Sobereva的最优先解答。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(十分强大、极为流行的量子化学波函数分析程序)
Google Scholar:https://scholar.google.com/citations?user=tiKE0qkAAAAJ
ResearchGate:https://www.researchgate.net/profile/Tian_Lu

6万

帖子

99

威望

5万

eV
积分
120134

管理员

公社社长

36#
 楼主 Author| 发表于 Post on 2025-5-3 00:01:31 | 只看该作者 Only view this author
zhouyulin 发表于 2025-5-2 23:50
老师,我想问一下优化结构对比结果后的RMSD怎么保存下来
如果你是指保存RMSD数值,复制粘贴就完了
如果你是指保存align后的结构,file - save coordinate
北京科音自然科学研究中心http://www.keinsci.com)致力于计算化学的发展和传播,长期开办极高质量的各种计算化学类培训:初级量子化学培训班中级量子化学培训班高级量子化学培训班量子化学波函数分析与Multiwfn程序培训班分子动力学与GROMACS培训班CP2K第一性原理计算培训班,内容介绍以及往届资料购买请点击相应链接查看。这些培训是计算化学从零快速入门以及进一步全面系统性提升研究水平的高速路!培训各种常见问题见《北京科音办的培训班FAQ》
欢迎加入北京科音微信公众号获取北京科音培训的最新消息,并避免错过网上有价值的计算化学文章!
欢迎加入人气极高、专业性特别强的理论与计算化学综合交流群思想家公社QQ群(群号见此链接),合计达一万多人。北京科音培训班的学员在群中可申请VIP头衔,提问将得到群主Sobereva的最优先解答。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(十分强大、极为流行的量子化学波函数分析程序)
Google Scholar:https://scholar.google.com/citations?user=tiKE0qkAAAAJ
ResearchGate:https://www.researchgate.net/profile/Tian_Lu

68

帖子

0

威望

455

eV
积分
523

Level 4 (黑子)

35#
发表于 Post on 2025-5-2 23:50:05 | 只看该作者 Only view this author
老师,我想问一下优化结构对比结果后的RMSD怎么保存下来

890

帖子

3

威望

1681

eV
积分
2631

Level 5 (御坂)

傻傻的木瓜

34#
发表于 Post on 2025-4-8 08:14:22 | 只看该作者 Only view this author
本帖最后由 Uus/pMeC6H4-/キ 于 2025-4-8 08:23 编辑

昨天我发的这帖就用此文的方法,简单展示了某个质子转移反应的反应物、过渡态、产物的叠合情况。

做些笔记:在VMD中启用File - Log Tcl Commands to Console后,可见从Extensions - Analysis打开RMSD Calculator窗口对应于命令行menu rmsd on指令,打开RMSD Trajectory Tool窗口对应于命令行menu rmsdtt on指令,打开RMSD Visualizer Tool窗口对应于命令行menu rmsdvt on指令。在窗口内操作时并不会在命令行显示对应指令,不过定义几个工具的原始tcl脚本可以在VMD目录下plugins/noarch/tcl/rmsd1.0或rmsdtt3.0或rmsdvt1.0文件夹内找到。实际上VMD自带的doc/ug.pdf用户指南的12.4节专门介绍了RMSD相关的tcl脚本编写技巧,核心指令就是measure fit $sel1 $sel2返回4x4矩阵,使得选区$sel1的原子改变坐标后与$sel2间的RMSD最小,以及measure rmsd $sel1 $sel2计算RMSD本身。measure fit有一个选项order可以指定作为参比的$sel2里的原子index列表。VMD的RMSD算法来自Kabsch (Acta Cryst. (1978) A34, 827-828)。

有一点不足之处是,这个方法得到的是视角(具体说就是http://bbs.keinsci.com/thread-17814-1-1.html涉及的若干矩阵)不变的前提下原子坐标的变换矩阵,但实际上有时可能需要的是体系的视角变换矩阵。比如对两个结构相似但朝向不同的分子用波函数信息计算实空间函数后将cube文件载入VMD,希望将原子位置和等值面同时在视觉上叠合,这就不能只做原子坐标的矩阵变换而不管格点数据了。不知道除了在产生波函数信息前就把分子对齐以外,有没有后处理上的解决方法。
√546=23.36664289109

1043

帖子

0

威望

4174

eV
积分
5217

Level 6 (一方通行)

33#
发表于 Post on 2023-10-4 16:54:33 来自手机 | 只看该作者 Only view this author
sobereva 发表于 2023-10-4 14:31
这是一个办法,也可以把距离值按从大到小排序成为矢量,比较两个矢量的差异,molclus的isostat比较结构相 ...

明白啦,谢谢!

6万

帖子

99

威望

5万

eV
积分
120134

管理员

公社社长

32#
 楼主 Author| 发表于 Post on 2023-10-4 14:31:34 | 只看该作者 Only view this author
granvia 发表于 2023-10-4 12:22
是不是比较距离矩阵的特征值就够了?

这是一个办法,也可以把距离值按从大到小排序成为矢量,比较两个矢量的差异,molclus的isostat比较结构相似性的时候就是这么做的
北京科音自然科学研究中心http://www.keinsci.com)致力于计算化学的发展和传播,长期开办极高质量的各种计算化学类培训:初级量子化学培训班中级量子化学培训班高级量子化学培训班量子化学波函数分析与Multiwfn程序培训班分子动力学与GROMACS培训班CP2K第一性原理计算培训班,内容介绍以及往届资料购买请点击相应链接查看。这些培训是计算化学从零快速入门以及进一步全面系统性提升研究水平的高速路!培训各种常见问题见《北京科音办的培训班FAQ》
欢迎加入北京科音微信公众号获取北京科音培训的最新消息,并避免错过网上有价值的计算化学文章!
欢迎加入人气极高、专业性特别强的理论与计算化学综合交流群思想家公社QQ群(群号见此链接),合计达一万多人。北京科音培训班的学员在群中可申请VIP头衔,提问将得到群主Sobereva的最优先解答。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(十分强大、极为流行的量子化学波函数分析程序)
Google Scholar:https://scholar.google.com/citations?user=tiKE0qkAAAAJ
ResearchGate:https://www.researchgate.net/profile/Tian_Lu

1043

帖子

0

威望

4174

eV
积分
5217

Level 6 (一方通行)

31#
发表于 Post on 2023-10-4 12:22:47 | 只看该作者 Only view this author
sobereva 发表于 2023-10-3 14:12
可以通过距离矩阵对比差异
molclus的isostat比较结构差异的时候就是基于距离矩阵对比的

是不是比较距离矩阵的特征值就够了?

6万

帖子

99

威望

5万

eV
积分
120134

管理员

公社社长

30#
 楼主 Author| 发表于 Post on 2023-10-3 14:12:29 | 只看该作者 Only view this author
granvia 发表于 2023-10-2 10:04
求教一下:如果两个结构的原子排列顺序不一致,有没有不手动而完全自动化计算RMSD的算法?

可以通过距离矩阵对比差异
molclus的isostat比较结构差异的时候就是基于距离矩阵对比的
北京科音自然科学研究中心http://www.keinsci.com)致力于计算化学的发展和传播,长期开办极高质量的各种计算化学类培训:初级量子化学培训班中级量子化学培训班高级量子化学培训班量子化学波函数分析与Multiwfn程序培训班分子动力学与GROMACS培训班CP2K第一性原理计算培训班,内容介绍以及往届资料购买请点击相应链接查看。这些培训是计算化学从零快速入门以及进一步全面系统性提升研究水平的高速路!培训各种常见问题见《北京科音办的培训班FAQ》
欢迎加入北京科音微信公众号获取北京科音培训的最新消息,并避免错过网上有价值的计算化学文章!
欢迎加入人气极高、专业性特别强的理论与计算化学综合交流群思想家公社QQ群(群号见此链接),合计达一万多人。北京科音培训班的学员在群中可申请VIP头衔,提问将得到群主Sobereva的最优先解答。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(十分强大、极为流行的量子化学波函数分析程序)
Google Scholar:https://scholar.google.com/citations?user=tiKE0qkAAAAJ
ResearchGate:https://www.researchgate.net/profile/Tian_Lu

1043

帖子

0

威望

4174

eV
积分
5217

Level 6 (一方通行)

29#
发表于 Post on 2023-10-2 10:04:17 | 只看该作者 Only view this author
求教一下:如果两个结构的原子排列顺序不一致,有没有不手动而完全自动化计算RMSD的算法?

41

帖子

0

威望

317

eV
积分
358

Level 3 能力者

28#
发表于 Post on 2023-7-20 20:12:43 | 只看该作者 Only view this author
sobereva 发表于 2023-7-20 19:22
用isostat默认的结构偏离标准就够了,没必要结合VMD

谢谢sob老师的指点

6万

帖子

99

威望

5万

eV
积分
120134

管理员

公社社长

27#
 楼主 Author| 发表于 Post on 2023-7-20 19:22:32 | 只看该作者 Only view this author
qcn1211 发表于 2023-7-20 17:32
请问sob老师,用molclus调动高斯计算筛选低能构象,结合VMD align, RMSD值为多少时说明结构相近?进而需要i ...

用isostat默认的结构偏离标准就够了,没必要结合VMD
北京科音自然科学研究中心http://www.keinsci.com)致力于计算化学的发展和传播,长期开办极高质量的各种计算化学类培训:初级量子化学培训班中级量子化学培训班高级量子化学培训班量子化学波函数分析与Multiwfn程序培训班分子动力学与GROMACS培训班CP2K第一性原理计算培训班,内容介绍以及往届资料购买请点击相应链接查看。这些培训是计算化学从零快速入门以及进一步全面系统性提升研究水平的高速路!培训各种常见问题见《北京科音办的培训班FAQ》
欢迎加入北京科音微信公众号获取北京科音培训的最新消息,并避免错过网上有价值的计算化学文章!
欢迎加入人气极高、专业性特别强的理论与计算化学综合交流群思想家公社QQ群(群号见此链接),合计达一万多人。北京科音培训班的学员在群中可申请VIP头衔,提问将得到群主Sobereva的最优先解答。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(十分强大、极为流行的量子化学波函数分析程序)
Google Scholar:https://scholar.google.com/citations?user=tiKE0qkAAAAJ
ResearchGate:https://www.researchgate.net/profile/Tian_Lu

41

帖子

0

威望

317

eV
积分
358

Level 3 能力者

26#
发表于 Post on 2023-7-20 17:32:33 | 只看该作者 Only view this author
请问sob老师,用molclus调动高斯计算筛选低能构象,结合VMD align, RMSD值为多少时说明结构相近?进而需要isostat重新归簇。

4104

帖子

4

威望

8869

eV
积分
13053

Level 6 (一方通行)

MOKIT开发者

25#
发表于 Post on 2023-5-29 10:46:12 | 只看该作者 Only view this author
qcn1211 发表于 2023-5-29 10:37
请问用VMD 计算出的RMSD值为多少才说明这两个结构相近?

说说你的体系多少个原子,你比较的是所有原子还是部分原子
自动做多参考态计算的程序MOKIT

6万

帖子

99

威望

5万

eV
积分
120134

管理员

公社社长

24#
 楼主 Author| 发表于 Post on 2023-5-29 10:40:00 | 只看该作者 Only view this author
qcn1211 发表于 2023-5-29 10:37
请问用VMD 计算出的RMSD值为多少才说明这两个结构相近?

没一刀切的标准。结合叠合的图像凭直觉判断最靠谱
北京科音自然科学研究中心http://www.keinsci.com)致力于计算化学的发展和传播,长期开办极高质量的各种计算化学类培训:初级量子化学培训班中级量子化学培训班高级量子化学培训班量子化学波函数分析与Multiwfn程序培训班分子动力学与GROMACS培训班CP2K第一性原理计算培训班,内容介绍以及往届资料购买请点击相应链接查看。这些培训是计算化学从零快速入门以及进一步全面系统性提升研究水平的高速路!培训各种常见问题见《北京科音办的培训班FAQ》
欢迎加入北京科音微信公众号获取北京科音培训的最新消息,并避免错过网上有价值的计算化学文章!
欢迎加入人气极高、专业性特别强的理论与计算化学综合交流群思想家公社QQ群(群号见此链接),合计达一万多人。北京科音培训班的学员在群中可申请VIP头衔,提问将得到群主Sobereva的最优先解答。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(十分强大、极为流行的量子化学波函数分析程序)
Google Scholar:https://scholar.google.com/citations?user=tiKE0qkAAAAJ
ResearchGate:https://www.researchgate.net/profile/Tian_Lu

41

帖子

0

威望

317

eV
积分
358

Level 3 能力者

23#
发表于 Post on 2023-5-29 10:37:53 | 只看该作者 Only view this author
请问用VMD 计算出的RMSD值为多少才说明这两个结构相近?

手机版 Mobile version|北京科音自然科学研究中心 Beijing Kein Research Center for Natural Sciences|京公网安备 11010502035419号|计算化学公社 — 北京科音旗下高水平计算化学交流论坛 ( 京ICP备14038949号-1 )|网站地图

GMT+8, 2025-8-15 06:14 , Processed in 0.198421 second(s), 25 queries , Gzip On.

快速回复 返回顶部 返回列表 Return to list