计算化学公社

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

[程序/脚本开发] 使用MULE在势能面/自由能面上寻找最低反应路径

[复制链接 Copy URL]

1149

帖子

6

威望

6629

eV
积分
7898

Level 6 (一方通行)

本帖最后由 fhh2626 于 2020-6-3 22:34 编辑

MUltidimensional Least Energy pathway (MULE)的原理其实就是Dijkstra算法,属于贪心算法的一种,具体数学内容就不多扯了。和其他软件相比,MULE有几个优点:
• 得到的一定是理论严谨的最低自由能路径,而其他方法,比如LFEP和String会陷入局部极小值
• 可以处理周期性反应坐标(比如二面角)
• 可以处理任意维度,现在网上能找到的大多数代码都只能处理二维路径
• 可以施加约束,从而实现A*算法,或者寻找第2/3/N条最低自由能路径,亦或者寻找经过某个点的最低自由能路径


MULE的windows主程序/源代码/示例文件都在压缩包里面提供,linux和mac用户需要自己编译mule.cpp文件


使用MULE需要写一个config文件,然后执行
./mule.exe config
软件就会自动生成一个.traj和一个.energy文件,分别对应路径和沿路径的能量/自由能变化


/example1文件夹中给出了config文件的写法
[mule]
directory       =   nanma_ref.pmf
initial         =    -156, 160
end             =      78, -58
pbc             =       1, 1
软件会读取执行目录的nanma_ref.pmf文件,因为这个文件是NAMD生成的,软件会自己读取里面记录的上下界和格点精度等信息,如果给出的文件是普通格式的话需要人为在config文件中提供这些信息。剩下的三行分别表示路径从(-156, 160)开始,到(78, -58)为止,反应坐标的两个维度都是周期性的,如果非周期性的反应坐标则将pbc设为0。找出来的路径如图中左上所示。


/example2文件夹中给出了找这个文件中第2低的反应路径的例子,config文件内容如下:
[mule]
directory       =   nanma_ref.pmf
initial         =   -156, 160
end             =     78, -58
lowerboundary   =   -180, -180
upperboundary   =    178, 178
width           =      2, 2
pbc             =     1, 1
writeExploredPoints   =   1
target          =     136, -32, -0.1, -0.1
这个例子中我提供了一个普通dat格式的文件,所以需要手动设置上下界和格点精度(lowerboundary, upperboundary and width)writeExploredPoints选项设为1会让程序输出探索过的所有点,一般来说设计A*算法时有些用,这里没用,不去管它。值得解释的是最后一行,这一行施加了一个Manhattan约束,具体数学公式不多说了。这里代表采用(0.1, 0.1)的力常数将探索过程推离(136, -32),其中(136, -32)是example1中得到的过渡态。这个约束会使得算法找到一条较为远离上一步中过渡态的路径,如右上子图所示


如果我们看这个自由能图,第一反应的话一定是找到从中间走的这个路径,example3给出了寻找这条路径的方法
[mule]
directory       =   nanma_ref.pmf
initial         =   -156, 160
end             =     78, -58
lowerboundary   =   -180, -180
upperboundary   =    178, 178
width           =      2, 2
pbc             =     1, 1
writeExploredPoints   =   1
target          =     0, 100, 0.1, 0.1
其他参数都一样,最后一行代表采用(0.1, 0.1)的力常数将探索过程拉向(0, 100)这个点(接近于这条路径的过渡态),得到的路径为左下图

右下的图对比了三条路径的自由能变化,可以看见第一条路径明显最有利





Reference:        Fu H., Chen H, Wang X, et al. Finding an Optimal Pathway on a Multidimensional Free-Energy Landscape. J. Chem. Inf. Model. 2020. doi:10.1021/acs.jcim.0c00279

MULE.zip (1.68 MB, 下载次数 Times of downloads: 303)







1.jpg (428.77 KB, 下载次数 Times of downloads: 187)

1.jpg

评分 Rate

参与人数
Participants 10
威望 +1 eV +33 收起 理由
Reason
ysss + 2 谢谢分享
iota + 4 牛!
hebrewsnabla + 3 好物!
snljty2 + 5 高级!
wzkchem5 + 5
丁越 + 5 好物!
壹零壹室掃地僧 + 2 赞!
sobereva + 1
ene + 5
minishotgun + 2 赞!

查看全部评分 View all ratings

2

帖子

0

威望

25

eV
积分
27

Level 2 能力者

28#
发表于 Post on 2025-7-3 09:37:27 | 只看该作者 Only view this author
fhh2626 发表于 2025-6-27 20:06
这是你提供的自由能面的精度决定的,你可以做一个插值或者平滑

非常感谢

1149

帖子

6

威望

6629

eV
积分
7898

Level 6 (一方通行)

27#
 楼主 Author| 发表于 Post on 2025-6-27 20:06:27 | 只看该作者 Only view this author
zhonghaoxuan 发表于 2025-6-23 10:48
老师您好,跟您请教一下,使用MULE的时候,我生成的反应路径上只有70个数据点,能量也只有70个,请问如何设 ...

这是你提供的自由能面的精度决定的,你可以做一个插值或者平滑

2

帖子

0

威望

25

eV
积分
27

Level 2 能力者

26#
发表于 Post on 2025-6-23 10:48:52 | 只看该作者 Only view this author
老师您好,跟您请教一下,使用MULE的时候,我生成的反应路径上只有70个数据点,能量也只有70个,请问如何设置CONFIG文件可以使数据点更多呢

1149

帖子

6

威望

6629

eV
积分
7898

Level 6 (一方通行)

25#
 楼主 Author| 发表于 Post on 2024-4-18 11:54:02 | 只看该作者 Only view this author
gehan 发表于 2024-4-18 10:48
再次感谢,那请问是怎么根据势能面上极值点的信息,定位到MD轨迹文件中的哪一帧结构呢,是根据极值点的对 ...

你知道极值点就知道CVs了,就可以去轨迹里筛选结构

30

帖子

0

威望

271

eV
积分
301

Level 3 能力者

24#
发表于 Post on 2024-4-18 10:48:20 | 只看该作者 Only view this author
fhh2626 发表于 2024-4-18 10:15
横坐标就是0-1自己插值

MD会输出轨迹文件啊,在轨迹文件里面找

再次感谢,那请问是怎么根据势能面上极值点的信息,定位到MD轨迹文件中的哪一帧结构呢,是根据极值点的对应的CV1和CV2,然后去fes.dat,或者该文章中的.pmf文件中找到对应的值,然后在该文件中找到对应的帧数吗。还是说是自己看MD轨迹动画去找到想要极值点结构呢。

1149

帖子

6

威望

6629

eV
积分
7898

Level 6 (一方通行)

23#
 楼主 Author| 发表于 Post on 2024-4-18 10:15:02 | 只看该作者 Only view this author
gehan 发表于 2024-4-18 10:01
谢谢您的回答,但是我使用MULE生成了.energy文件里是只有一列能量的数值,没有横坐标哎,不知道怎么设置 ...

横坐标就是0-1自己插值

MD会输出轨迹文件啊,在轨迹文件里面找

30

帖子

0

威望

271

eV
积分
301

Level 3 能力者

22#
发表于 Post on 2024-4-18 10:01:24 | 只看该作者 Only view this author
fhh2626 发表于 2024-4-16 14:22
origin可以添加图层,添加一个图层画散点图就行了

右下角的图是以路径为坐标,如果你用的MULE的话会直 ...

谢谢您的回答,但是我使用MULE生成了.energy文件里是只有一列能量的数值,没有横坐标哎,不知道怎么设置路径为坐标。还有一个问题,请问有没有什么办法,找到势能面中的极小值点以及过渡态的位置后,怎么才能找到对应的结构呢?

1149

帖子

6

威望

6629

eV
积分
7898

Level 6 (一方通行)

21#
 楼主 Author| 发表于 Post on 2024-4-16 14:22:36 | 只看该作者 Only view this author
gehan 发表于 2024-4-16 10:12
请问楼主,左上,左下和右上三个图是怎么绘制出的,二维等值面可以使用origin生成,图中的黑色路径的线是怎 ...

origin可以添加图层,添加一个图层画散点图就行了

右下角的图是以路径为坐标,如果你用的MULE的话会直接输出这个文件

30

帖子

0

威望

271

eV
积分
301

Level 3 能力者

20#
发表于 Post on 2024-4-16 10:12:19 | 只看该作者 Only view this author
请问楼主,左上,左下和右上三个图是怎么绘制出的,二维等值面可以使用origin生成,图中的黑色路径的线是怎么在图中显示出来的,还有就右下的自由能变化曲线的横坐标是取什么为变量的呢?谢谢

1149

帖子

6

威望

6629

eV
积分
7898

Level 6 (一方通行)

19#
 楼主 Author| 发表于 Post on 2024-1-28 11:44:41 | 只看该作者 Only view this author
mz121 发表于 2024-1-26 13:35
这个自由能面是通过构建马尔科夫链得到的,现在我已经填满了空间,但mule查找的自由能路径上的自由能值大 ...

你的dat/pmf文件是等间隔的吗,你看看例子中的文件,弄成一样的格式应该就行

3

帖子

0

威望

21

eV
积分
24

Level 1 能力者

18#
发表于 Post on 2024-1-26 13:35:04 | 只看该作者 Only view this author
本帖最后由 mz121 于 2024-1-26 13:36 编辑
fhh2626 发表于 2024-1-26 09:36
你写个脚本,让自由能面填满整个空间,没有值的地方都填个很大的数就可以(Colvars和Plumed输出的自由能 ...

这个自由能面是通过构建马尔科夫链得到的,现在我已经填满了空间,但mule查找的自由能路径上的自由能值大多数为0,我自由能面中没有0值,不知道mule为何会查找出0值?

微信截图_20240126133519.png (281.46 KB, 下载次数 Times of downloads: 144)

微信截图_20240126133519.png

3

帖子

0

威望

21

eV
积分
24

Level 1 能力者

17#
发表于 Post on 2024-1-26 13:27:13 | 只看该作者 Only view this author
fhh2626 发表于 2024-1-26 09:36
你写个脚本,让自由能面填满整个空间,没有值的地方都填个很大的数就可以(Colvars和Plumed输出的自由能 ...

这个自由能是通过马尔可夫链得出的。我写了个脚本,填充了空间,但mule查找出来的自由能路径上的自由能值为0。我自由能文件中并未有0值,不知道为何mule会得到0值?

1149

帖子

6

威望

6629

eV
积分
7898

Level 6 (一方通行)

16#
 楼主 Author| 发表于 Post on 2024-1-26 09:39:25 | 只看该作者 Only view this author
elpa 发表于 2023-2-16 19:58
最后几列具体是哪几列?应该一共5列,最后两列是微分。是不是只留前三列?能不能给我一个参考示例,fes.d ...

PMF和dat是一个格式,就是

CV1 CV2 ... CVn deltaG

评分 Rate

参与人数
Participants 1
eV +3 收起 理由
Reason
elpa + 3 谢谢

查看全部评分 View all ratings

1149

帖子

6

威望

6629

eV
积分
7898

Level 6 (一方通行)

15#
 楼主 Author| 发表于 Post on 2024-1-26 09:38:48 | 只看该作者 Only view this author
youyno 发表于 2023-2-16 16:46
把最后几列不需要的删掉后再改个后缀就行,关键是怎么确定极小值点

可以写个简单的贪心算法程序找,另外即使不选极小值点,MULE找出来的路径肯定也会经过附近的极小值点,直接从输出文件中找就行

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

GMT+8, 2025-8-14 07:03 , Processed in 0.206345 second(s), 31 queries , Gzip On.

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