计算化学公社
标题: 谈谈记录化学体系结构的mol2文件 [打印本页]
作者Author: sobereva 时间: 2022-12-31 22:03
标题: 谈谈记录化学体系结构的mol2文件
谈谈记录化学体系结构的mol2文件
Introduction to the mol2 format for recording chemical structures
文/Sobereva@北京科音 2022-Dec-31
1 mol2文件简介
mol2文件是计算化学领域非常常用的记录分子结构的格式,被很多程序所支持和利用。例如VMD、GaussView、Chem3D、OpenBabel、AmberTools里的Antechamber等程序都可以导出mol2文件,笔者开发的Multiwfn(http://sobereva.com/multiwfn)可以基于mol2文件计算EEM原子电荷,笔者开发的Sobtop(http://sobereva.com/soft/Sobtop)可以基于mol2文件产生GROMACS的拓扑文件,等等。
相对于非常常见的在《谈谈记录化学体系结构的xyz文件》(http://sobereva.com/477)中介绍的xyz格式,mol2格式关键优点在于可以记录成键信息,即谁与谁成什么形式的键,这对于判断原子所处的化学环境非常重要,比如Sobtop需要有这样的成键信息才能自动指认GAFF原子类型,Multiwfn计算EEM原子电荷时需要有这样的信息才能判断各个原子要用的EEM参数。另外mol2文件还可以记录其它诸多信息(距离约束、可旋转的键、原子类型和电荷等等),对于分子模拟、QSAR、化学信息学等一些方面有特殊的意义。
mol2格式略微复杂,不同程序产生的mol2文件有所出入,有的程序产生的mol2文件还不规矩,导致经常有人由于用的mol2文件有问题而无法用Sobtop和Multiwfn等程序正常处理,甚至导致程序崩溃。我遂觉得有必要写一篇文章介绍一下mol2格式,便于读者了解怎么读取mol2文件的信息、判断自己手头的mol2文件是否规范,以及拿到不标准的mol2文件时怎么修改。
mol2文件是文本格式,包含大量的字段,每个字段各有各的用处和定义规范。mol2文件的详细说明可以下载此文档查看:http://sobereva.com/attach/655/Tripos_Mol2_File_Format.pdf。这些字段并不是全都需要出现的,常见的字段只有几个而已,每个字段涉及的信息中往往也只有其中少部分会经常涉及到,将在下文进行介绍,若想更全面详细了解mol2格式可阅读上述文档,共54页。
本文提到的程序的版本是VMD 1.9.3、Multiwfn 3.8(dev, 2022-Dec-18)、Sobtop 1.0(dev3.1)、GaussView 6.0.16、OpenBabel 3.1.1。其它版本可能与本文相同也可能不同。
2 mol2文件的例子和解读
OpenBabel程序产生的mol2格式相对来说是属于比较规矩的,这里结合OpenBabel程序产生的苯甲醛的mol2文件的内容进行讲解。我先用一个可视化程序画了一个苯甲醛的结构,保存为了benzaldehyde.pdb文件,然后用obabel benzaldehyde.pdb -O benzaldehyde.mol2命令用OpenBabel把pdb格式转成了mol2格式的benzaldehyde.mol2文件,其内容如下,GaussView载入此文件后显示的结构图也给出了以便于对照
(, 下载次数 Times of downloads: 57)
@<TRIPOS>MOLECULE
benzaldehyde.pdb
14 14 0 0 0
SMALL
GASTEIGER
@<TRIPOS>ATOM
1 C 1.9920 0.4700 -0.0000 C.2 0 UNK0 0.1502
2 C 0.5340 0.2150 -0.0000 C.ar 0 UNK0 0.0142
3 C -0.3610 1.2920 -0.0000 C.ar 0 UNK0 -0.0515
4 C -1.7360 1.0600 0.0000 C.ar 0 UNK0 -0.0611
5 C -2.2160 -0.2520 0.0000 C.ar 0 UNK0 -0.0617
6 C -1.3250 -1.3310 -0.0000 C.ar 0 UNK0 -0.0611
7 C 0.0460 -1.1010 -0.0000 C.ar 0 UNK0 -0.0515
8 O 2.8450 -0.3960 0.0000 O.2 0 UNK0 -0.2957
9 H 2.2730 1.5470 0.0000 H 0 UNK0 0.1081
10 H 0.0250 2.3090 -0.0000 H 0 UNK0 0.0624
11 H -2.4300 1.8950 0.0000 H 0 UNK0 0.0618
12 H -3.2860 -0.4350 0.0000 H 0 UNK0 0.0618
13 H -1.7060 -2.3480 -0.0000 H 0 UNK0 0.0618
14 H 0.7610 -1.9170 -0.0000 H 0 UNK0 0.0624
@<TRIPOS>BOND
1 1 2 1
2 1 8 2
3 1 9 1
4 2 3 ar
5 2 7 ar
6 3 4 ar
7 3 10 1
8 4 5 ar
9 4 11 1
10 5 6 ar
11 5 12 1
12 6 7 ar
13 6 13 1
14 7 14 1
首先要知道mol2文件里以#作为第一列的是注释行,空行也被完全无视。mol2文件是自由格式,因此空格数目完全随意。
第一列以@开头的叫做字段,从上面的benzaldehyde.mol2可见,当前文件有@<TRIPOS>MOLECULE、@<TRIPOS>ATOM、@<TRIPOS>BOND三个字段。一般来说这三个字段都是必须出现的,一起提供了描述一个分子最起码的信息。
2.1 MOLECULE字段
@<TRIPOS>MOLECULE字段记录了体系的基本信息,包括:
第1行:体系的名字。可见OpenBabel把转换出mol2文件的源文件的名字benzaldehyde.pdb当做了当前体系的名字
第2行:五个数字分别是体系中的原子数、化学键数、子结构数、特征数、set数。对于单纯记录体系结构信息,只要提供前两者就够了,后三个可以省略。所谓的子结构是指体系中的一个部分,比如每个分子、每个残基、蛋白质的每条链等等都可以在@<TRIPOS>SUBSTRUCTURE字段里定义为一个子结构。所谓的set是指基于体系中的一些原子/键/子结构根据特定规则和需要定义的集合,可以在@<TRIPOS>SET里具体定义。
第3行:体系的类型。可以为SMALL(小分子)、BIOPOLYMER、PROTEIN、NUCLEIC_ACID、SACCHARIDE
第4行:原子电荷。如果mol2文件没记录原子电荷信息这里就为NO_CHARGES。而在产生当前benzaldehyde.mol2文件的时候OpenBabel自动计算了Gasteiger电荷,因此此处写的是GASTEIGER。还可以为MULLIKEN_CHARGES(Mulliken电荷)、MMFF94_CHARGES(MMFF94力场定义的电荷)等等,不同种类电荷都有固定名字。如果记录的原子电荷是比如Multiwfn算的ADCH、RESP、1.2*CM5等电荷,在mol2格式规范中没有对应的名字,则这里应当写USER_CHARGES。
2.2 ATOM字段
@<TRIPOS>ATOM字段每一行定义一个原子的信息,每一列记录的信息为:
(1)原子序号(整数)
(2)原子名(字符串)
(3)X坐标(埃)
(4)Y坐标(埃)
(5)Z坐标(埃)
(6)原子类型(atom type。字符串)
(7)原子所属的子结构序号(整数),可省略
(8)原子所述的子结构名字(字符串),可省略
(9)原子电荷(浮点数),可省略
原子名部分可以为比如C2、H4等等,完全随意。记录生物分子结构时通常用IUPAC定义的各种残基中的原子名。
原子类型部分可以记录做分子模拟用的力场中此原子实际对应的原子类型。mol2格式自己也有一套原子类型定义,见前述的Tripos_Mol2_File_Format.pdf文档的末尾,比如sp3杂化的碳的原子类型是C.3,C.ar代表芳香碳,Any代表任意,Hal泛指卤素,Cl代表氯,Ca代表钙,H代表氢,H.spc特指SPC水模型的氢,LP代表孤对电子(lone pair),Du代表虚原子(dummy),等等。
一定要特别注意,mol2格式虽然定义了一大堆字段,但(居然)没有一个地方是专门用来记录元素的,这在我来看是mol2格式的严重不足!!!mol2记录的原子名和原子类型信息可以与元素名相同也可以不同,不同程序产生的mol2文件的情况各有不同。例如如此例可见,OpenBabel产生的这个mol2文件里原子名恰等于元素名,原子类型是根据mol2格式自己定义的原子类型指认的。而在GaussView产生的mol2文件中,原子名是给元素名后面加了数字(因此不会有重名的原子),而原子类型恰等于元素名。由于情况混乱,所以一个程序在读取mol2文件的时候并没有严格的办法能准确判断元素,只能靠猜。Multiwfn和Sobtop在读取mol2文件时是根据原子类型的字符串判断元素的:如果字符串中没有.,就直接将之当做元素名来判断元素;如果有.,比如是C.3,就把.前面的内容当做元素名判断元素。因此,读者应该知道Multiwfn和Sobtop没判断对元素时该怎么办了,最简单的做法就是手动修改mol2文件以让@<TRIPOS>ATOM字段每一行的第6列对应元素名。
由于子结构信息在原子电荷前面,因此即便你不想定义原子所属的子结构信息而只想定义原子电荷,也必须随便写上子结构序号和子结构名字来占位,比如此例用0 UNK0来占位。
2.3 BOND字段
@<TRIPOS>BOND字段每一行定义一个键的信息,其每一列记录的信息为:
(1)键的序号(整数)
(2)第1个原子的序号
(3)第2个原子的序号
(4)键的类型
键的类型有以下这些
• 1 = 单键
• 2 = 双键
• 3 = 三键
• am = 酰胺的N-C键(这种键有一定pi共轭作用,这是为什么mol2格式里特意用am来与单键区分)
• ar = 芳香环(aromatic)上的键,以下简称芳香键
• du = 虚键
• un = 未知/无法判断
• nc = 不相连(俩原子不成键就没必要在BOND字段出现,但可以靠nc强调某两个原子间就是没成键)
绝大多数程序产生的mol2文件里没有du、un、nc。
有的程序产生的键的类型名不规矩。比如GaussView对于芳香环上的键(具体来说,是图形窗口里看到一个实线+一个虚线的那种键)用的类型是Ar,但按照mol2规范应当是ar,因此这会导致一些程序无法识别(而Sobtop在读取时已经考虑到了GaussView的这个bug,因此用户不用自己做替换)。GaussView还自行给mol2格式做了扩展,把图形窗口里看到是一个虚键的那种键记录为Wk代表Weak,而这可能导致很多程序无法正常识别和载入。
不同程序对成键的指认也往往有很大不同。比如甲酰胺,Avogadro产生的它的mol2里把N-C键记录为am,严格符合mol2格式的要求。而GaussView则无法将之记录为am,而是可能记录为单键也可能记录为Ar(取决于当前图形窗口里显示的成键方式)。另外,我之前在《谈谈原子间是否成键的判断问题》(http://sobereva.com/414)中说过GaussView是根据原子半径和原子间距离判断成键形式的,导致很有可能判断出的成键方式很不“经典”,甚至很违背化学常识,比如可能显示一个碳原子连着一个双键和两个单+虚键。而如果把结构保存成pdb然后再用OpenBabel转成mol2格式,则成键方式就很经典了,因为OpenBabel能够自动让成键方式满足经典Lewis式且同时识别芳香区域。
VMD程序里如果载入了xyz、gro、pdb等不含成键关系信息的文件(虽然pdb有CONECT字段,但VMD不利用),保存出的mol2文件将没有BOND字段,明显是不符合规范的。在同时载入拓扑文件以提供拓扑信息后,保存出的mol2文件才是有BOND字段的(与此同时也有了原子类型、原子电荷、原子所属的残基号和残基名信息)。VMD如果载入的是mol2文件,也可以从中获取成键信息,使得保存出的mol2文件里也有BOND部分。但是VMD并不能记录芳香键,而且它自己也没有像OpenBabel那样根据几何结构和元素就能判断出芳香区域的能力,因此即便载入的是本例的benzaldehyde.mol2,保存出的mol2里苯环上的键也都会简单地记录成单键。
还值得顺带一提的是有个叫mol的格式,介绍见https://en.wikipedia.org/wiki/Chemical_table_file。不要将它和mol2混淆,二者格式截然不同。mol也能像mol2一样记录键的存在性及其类型。GaussView产生的mol文件中会把芳香键用4来记录,这正是mol标准格式中的芳香键的类型序号,而OpenBabel在产生mol文件时则会把芳香环描述成单双键交替的形式。
3 关于记录晶胞信息
mol2文件通常用来记录孤立体系,但实际上此格式也定义了记录晶胞信息的字段@<TRIPOS>CRYSIN,在其下一行写晶胞的a、b、c三个边长(埃)以及alpha、beta、gamma夹角(度),每个值之间以逗号分隔。例如:
@<TRIPOS>CRYSIN
3.785,3.785,9.514,90,90,90
Sobtop和Multiwfn在载入mol2文件时都会试图载入晶胞信息以考虑周期性。GaussView有很好用的对周期性体系建模的功能,但即便在图形窗口里能看到晶胞边框,代表当前是周期性体系,其保存出的mol2文件里也没有@<TRIPOS>CRYSIN字段,因此需要自行补充。在VMD 1.9.3里,即便当前有晶胞信息,保存的mol2文件里也没有以上字段。VMD 1.9.3载入mol2时也不会从以上字段中载入晶胞信息。
作者Author: 牧生 时间: 2023-1-1 09:00
(, 下载次数 Times of downloads: 58)
我一直以为gview导出的mol2格式最规范,原来不是。
“OpenBabel程序产生的mol2格式相对来说是属于比较规矩的”,那么有没有例外呢?可不可以认为在当前,OpenBabel就是最规范的?
作者Author: sobereva 时间: 2023-1-2 06:29
规范在本文和那个mol2文件介绍文档都介绍了,对照一下就知道规不规范。
OpenBabel产生的mol2文件在常见程序中属于较高的(但不敢说最高),比如原子类型部分都能按照mol2的规范识别并记录。但也不能说绝对完美(其实也和mol2格式定义得太细碎有关),比如酰胺的N-C键的类型就没严格按照mol2格式的规范记录成am,不过这不算什么大问题。
作者Author: lyj714 时间: 2023-1-2 11:33
本帖最后由 lyj714 于 2023-1-2 11:45 编辑
其实vmd(包括1.9.3 win32等各个版本)是可以利用pdb格式中的CONECT的,这部分原子连接作为成键搜索以后的补充。不过这个部分必须每一行只能定义一个键,然后每一个键都需要单独写一行CONECT,并且每一行的原子序号必须是递增的(比如1和2原子有键,就写成1 2而不是2 1),这样vmd就可以把CONECT部分作为键补充部分加到结构中显示。比如给一个水加上H-H的键显示:
- ATOM 1 OW SOL A 4 17.180 36.520 6.210 1.00 0.00 O
- ATOM 2 HW1 SOL A 4 16.220 36.560 6.250 1.00 0.00 H
- ATOM 3 HW2 SOL A 4 17.430 37.160 5.560 1.00 0.00 H
- CONECT 2 3
复制代码
(, 下载次数 Times of downloads: 59)
作者Author: sobereva 时间: 2023-1-3 03:47
利用pdb里的CONECT信息,用此文的做法比较方便
使VMD根据pdb文件中的CONECT字段设定原子连接关系
http://sobereva.com/121
作者Author: hanlan8702 时间: 2023-1-11 12:36
请问文中提到得“还可以为MULLIKEN_CHARGES(Mulliken电荷)、MMFF94_CHARGES(MMFF94力场定义的电荷)等等,不同种类电荷都有固定名字。”这部分用openbabel应该如何实现呢?谢谢。
作者Author: sobereva 时间: 2023-1-12 06:52
我不知道你说的实现是什么意思
如果是OpenBabel计算的话,它算不了Mulliken电荷。计算MMFF94电荷可以在Openbabel转换时加上--partialcharge MMFF94
作者Author: hanlan8702 时间: 2023-1-12 07:17
好的,谢谢老师。
作者Author: haijiayuan 时间: 2023-2-19 13:37
sob老师,我想把lammps的data文件转换成mol2格式(带成键信息),请问有什么推荐的方法吗?
作者Author: sobereva 时间: 2023-2-20 02:25
我不用lammps
作者Author: haijiayuan 时间: 2023-2-22 10:26
好的,谢谢
作者Author: 秋心 时间: 2024-7-22 09:35
sob老师您好,我在使用OpenBabel产生二甲胺的mol2文件时发现产生的mol2文件只有一种氢原子类型,很明显N上的氢原子与C上的氢原子化学环境不同,但是OpenBabel却归为了一类,看OpenBabel的mol2文件格式里面确实没法区分N和C上氢原子。我想问一下如果我就用mol2产生的原子类型使用NWchem去做ISA计算得到该分子原子多极矩,弥散系数等,最后的结果可信吗?
作者Author: sobereva 时间: 2024-7-23 03:41
我不清楚你说的那是什么功能,我不评论。
如果要得到原子的色散系数,用下文的方法,根本不牵扯到原子类型
使用Multiwfn计算原子的C6色散系数
http://sobereva.com/709(http://bbs.keinsci.com/thread-45686-1-1.html)
用Multiwfn的主功能15基于模糊空间,或者基于主功能17基于AIM原子空间,去计算原子多极矩也完全不牵扯到原子类型。看Multiwfn手册4.15节和4.17节的相关例子。
作者Author: 秋心 时间: 2024-7-23 09:24
感谢sob老师。学习到了,原来multiwfn还可以计算色散系数啊,看来还是有必要再深入了解multiwfn的功能。
作者Author: Uus/pMeC6H4-/キ 时间: 2025-2-11 11:54
社长你好,我还有点疑惑,BOND字段里ar成键类型的适用范围有多广呢?比如,吲哚的六元环上的键全是ar,但独属于五元环的氮碳键、碳碳键应该定为ar么?[18]轮烯的碳碳键,孤立环丁二烯、[Fe(CO)3(C4H4)]配合物里环丁二烯的碳碳键,胍正离子的氮碳键呢?
虽然用OpenBabel产生mol2文件的确可以给上述问题一个明确答案(比如吲哚的五元环也全是ar),量子化学计算也不会关心这些形式键级,但是我感觉创建拓扑文件时需要关心,毕竟mol2文件的成键类型被程序读取后可能用来判断原子类型和成键参数。据说sobtop及其附带的atomtype同样不看形式键级只取连接关系,不知道其他的程序和力场如何。
作者Author: sobereva 时间: 2025-2-12 05:27
属于五元环的氮碳键、碳碳键也算ar,毕竟在pi共轭范围内。算不算ar没有严格的确切的界定标准
作者Author: Uus/pMeC6H4-/キ 时间: 2025-2-12 08:14
那还是看实际的体系和力场对形式键级的需要来定义成键类型好了。
这篇博文2.2节说到mol2的原子类型定义有“Du代表虚原子(dummy)”,但是GaussView能读取却不能保存含有虚原子的mol2文件,而且sobtop也没法识别。如果能支持虚原子的话,或许有助于定义一些特殊的“成键项”。比如类似http://sobereva.com/635的二茂铁拓扑文件,但俩茂环中心各有一虚原子Du,约束与各自环上的碳原子的Du-C“键长”,这样可以产生Du-Fe“键长”项、Du1-Fe-Du2“键角”项和C1-Du1-Du2-C2跨过铁的“二面角”项描述俩环的整体相对运动。
| 欢迎光临 计算化学公社 (http://bbs.keinsci.com/) |
Powered by Discuz! X3.3 |