增强采样软件PLUMED的安装与使用

技术背景

增强采样（Enhanced Sampling）是一种在分子动力学模拟中常用的技术，其作用是帮助我们更加快速的在时间轴上找到尽可能多的体系结构及其对应的能量。比如一个氢气的燃烧反应，在中间过程中会产生众多的反应产物，但是我们光从结果来看的话，就是从\(H_2\)加\(O_2\)，经过燃烧以后变成了\(H_2O\)，那么中间的过程就被我们全部忽略了。实际上化学反应是一个整体的过程，是非常复杂的，如果要控制一个化学反应的过程，或者制备某个中间态的产物，或者需要提升最终预期产物的占比，我们不得不考虑整个化学反应的路径。

由于分子动力学模拟是基于牛顿力学的，因此中间发生的位移和碰撞我们都认为是连续的过程，但是为了在时间轴上更快的读取数据，我们不能无限制的对时间进行分割，只能通过对时间进行采样来获得一个近似的过程。最典型的就是均匀采样，比如每隔10ps采一个样本点。但是这种方法在保障计算效率的前提下，很容易忽略了中间过程中出现时间极其短暂的一个体系状态。因此就需要使用到增强抽样的方法，对于高简并度的状态，我们降低其被采样的概率，而对于低简并度的状态，我们提升其被采到的概率。常见的方法有：Meta Dynamics、VES和ITS等。这里我们探索一下分子动力学模拟软件PLUMED的安装，该软件已经集成了很大一部分的CV和增强采样的方法。

PLUMED的下载与安装

首先访问其官方下载网站找到一个合适的版本进行下载，比如这里我是直接下载的最新的版本：

下载完成后可以用tar -xvf plumed-2.7.1.tgz指令来进行解压缩，可以看到解压后的目录如下所示：

(base) dechin@ubuntu2004:~/projects/gitlab/dechin/src/plumed$ ll
总用量 103736
drwxrwxr-x  3 dechin dechin      4096 7月  12 09:21 ./
drwxrwxr-x 10 dechin dechin      4096 7月  12 09:20 ../
drwxrwxr-x 18 dechin dechin      4096 4月  16 14:48 plumed-2.7.1/
-rw-rw-r--  1 dechin dechin 106210796 7月  12 09:21 plumed-2.7.1.tgz
(base) dechin@ubuntu2004:~/projects/gitlab/dechin/src/plumed$ cd plumed-2.7.1/
(base) dechin@ubuntu2004:~/projects/gitlab/dechin/src/plumed/plumed-2.7.1$ ll
总用量 504
drwxrwxr-x 18 dechin dechin   4096 4月  16 14:48 ./
drwxrwxr-x  3 dechin dechin   4096 7月  12 09:21 ../
drwxrwxr-x  6 dechin dechin   4096 4月  16 14:48 astyle/
-rw-rw-r--  1 dechin dechin     71 4月  16 14:48 .astyle.options
drwxrwxr-x  2 dechin dechin   4096 4月  16 14:48 CHANGES/
-rw-rw-r--  1 dechin dechin     94 4月  16 14:48 .codecov.yml
drwxrwxr-x  4 dechin dechin   4096 4月  16 14:48 conda/
-rwxrwxr-x  1 dechin dechin 323087 4月  16 14:48 configure*
-rw-rw-r--  1 dechin dechin  43692 4月  16 14:48 configure.ac
-rw-rw-r--  1 dechin dechin   7639 4月  16 14:48 COPYING.LESSER
drwxrwxr-x  2 dechin dechin   4096 4月  16 14:48 developer-doc/
drwxrwxr-x  2 dechin dechin   4096 4月  16 14:48 docker/
drwxrwxr-x  3 dechin dechin   4096 4月  16 14:48 .github/
-rw-rw-r--  1 dechin dechin    179 4月  16 14:48 .gitignore
-rw-rw-r--  1 dechin dechin    245 4月  16 14:48 .lgtm.yml
drwxrwxr-x  2 dechin dechin   4096 4月  16 14:48 macports/
-rw-rw-r--  1 dechin dechin   2485 4月  16 14:48 Makefile
-rw-rw-r--  1 dechin dechin   1200 4月  16 14:48 Makefile.conf.in
drwxrwxr-x  7 dechin dechin   4096 4月  16 14:48 patches/
-rw-rw-r--  1 dechin dechin    527 4月  16 14:48 PEOPLE
drwxrwxr-x  3 dechin dechin   4096 4月  16 14:48 python/
-rw-rw-r--  1 dechin dechin   8146 4月  16 14:48 README.md
drwxrwxr-x 26 dechin dechin   4096 4月  16 14:48 regtest/
-rwxrwxr-x  1 dechin dechin   5397 4月  16 14:48 release.sh*
drwxrwxr-x  2 dechin dechin   4096 4月  16 14:48 scripts/
-rw-rw-r--  1 dechin dechin    319 4月  16 14:48 sourceme.sh.in
drwxrwxr-x 45 dechin dechin   4096 4月  16 14:48 src/
-rw-rw-r--  1 dechin dechin      0 4月  16 14:48 stamp-h.in
drwxrwxr-x  6 dechin dechin   4096 4月  16 14:48 test/
drwxrwxr-x  2 dechin dechin   4096 4月  16 14:48 .travis/
-rw-rw-r--  1 dechin dechin  10558 4月  16 14:48 .travis.yml
drwxrwxr-x  4 dechin dechin   4096 4月  16 14:48 user-doc/
-rw-rw-r--  1 dechin dechin    322 4月  16 14:48 VERSION
drwxrwxr-x  2 dechin dechin   4096 4月  16 14:48 vim/

直接在这个plumed-2.7.1（注意对应自己安装的版本）目录下执行如下的指令进行安装：

$ ./configure --prefix=/usr/local
$ make -j 4
$ sudo make install

安装验证

安装完成后，可以在命令行中直接使用plumed指令，可以用如下帮助指令，确认软件是否被安装成功：

(base) dechin@ubuntu2004:~/projects/gitlab/dechin/src/plumed/examples$ plumed help
Usage: plumed [options] [command] [command options]
  plumed [command] -h|--help: to print help for a specific command
Options:
  [help|-h|--help]          : to print this help
  [--is-installed]          : fails if plumed is not installed
  [--has-mpi]               : fails if plumed is running without MPI
  [--has-dlopen]            : fails if plumed is compiled without dlopen
  [--load LIB]              : loads a shared object (typically a plugin library)
  [--standalone-executable] : tells plumed not to look for commands implemented as scripts
Commands:
  plumed completion : dump a function usable for programmable completion
  plumed driver : analyze trajectories with plumed
  plumed driver-float : analyze trajectories with plumed (single precision version)
  plumed gen_example : construct an example for the manual that users can interact with
  plumed gentemplate : print out a template input for a particular action
  plumed info : provide informations about plumed
  plumed kt : print out the value of kT at a particular temperature
  plumed manual : print out a description of the keywords for an action in html
  plumed pathtools : print out a description of the keywords for an action in html
  plumed pdbrenumber : Modify atom numbers in a PDB, possibly using hybrid-36 coding
  plumed pesmd : Langevin dynamics on PLUMED energy landscape
  plumed simplemd : run lj code
  plumed sum_hills : sum the hills with  plumed
  plumed patch : patch an MD engine
  plumed vim2html : convert plumed input file to colored html using vim syntax
  plumed selector : create lists of serial atom numbers
  plumed config : inquire plumed about how it was configure
  plumed newcv : create a new collective variable from a template
  plumed mklib : compile a .cpp file into a shared library
  plumed partial_tempering : scale parameters in a gromacs topology to implement solute or partial tempering

简单示例

使用plumed，需要有一个完整的反应路径文件，这里我们可以用wget直接下载plumed仓库中的示例轨迹文件：

(base) dechin@ubuntu2004:~/projects/gitlab/dechin/src/plumed/examples$ wget https://plumed.github.io/doc-v2.4/user-doc/html/tutorial-resources/trieste-1.tar.gz
--2021-07-12 10:18:18--  https://plumed.github.io/doc-v2.4/user-doc/html/tutorial-resources/trieste-1.tar.gz
正在解析主机 plumed.github.io (plumed.github.io)... ::1, 185.199.109.153, 185.199.108.153, ...
正在连接 plumed.github.io (plumed.github.io)|::1|:443... 失败：拒绝连接。
正在连接 plumed.github.io (plumed.github.io)|185.199.109.153|:443... 已连接。
已发出 HTTP 请求，正在等待回应... 301 Moved Permanently
位置：https://www.plumed.org/doc-v2.4/user-doc/html/tutorial-resources/trieste-1.tar.gz [跟随至新的 URL]
--2021-07-12 10:18:19--  https://www.plumed.org/doc-v2.4/user-doc/html/tutorial-resources/trieste-1.tar.gz
正在解析主机 www.plumed.org (www.plumed.org)... ::1, 185.199.109.153, 185.199.108.153, ...
正在连接 www.plumed.org (www.plumed.org)|::1|:443... 失败：拒绝连接。
正在连接 www.plumed.org (www.plumed.org)|185.199.109.153|:443... 已连接。
已发出 HTTP 请求，正在等待回应... 200 OK
长度： 335374 (328K) [application/gzip]
正在保存至: “trieste-1.tar.gz”

trieste-1.tar.gz    100%[==================>] 327.51K  --.-KB/s    用时 0.04s 

2021-07-12 10:18:20 (8.02 MB/s) - 已保存 “trieste-1.tar.gz” [335374/335374])

同样的执行解压：

(base) dechin@ubuntu2004:~/projects/gitlab/dechin/src/plumed/examples$ ll
总用量 336
drwxrwxr-x 2 dechin dechin   4096 7月  12 10:18 ./
drwxrwxr-x 4 dechin dechin   4096 7月  12 09:38 ../
-rw-rw-r-- 1 dechin dechin 335374 11月 18  2020 trieste-1.tar.gz
(base) dechin@ubuntu2004:~/projects/gitlab/dechin/src/plumed/examples$ tar -xvf trieste-1.tar.gz 
trieste-1/
trieste-1/.solutions/
trieste-1/.solutions/ref-rna.pdb
trieste-1/.solutions/plumed-ex2c.dat
trieste-1/.solutions/plumed-ex2.dat
trieste-1/.solutions/plumed-ex2b.dat
trieste-1/.solutions/plumed-ex1b.dat
trieste-1/.solutions/plumed-ex1.dat
trieste-1/traj-broken.xtc
trieste-1/ref.pdb
trieste-1/traj-whole.xtc

解压后可以在目录下看到一个pdb的配置文件和两个xtc的轨迹文件，这个是用gromacs生成的数据文件格式：

(base) dechin@ubuntu2004:~/projects/gitlab/dechin/src/plumed/examples/trieste-1$ ll
总用量 752
drwxr-xr-x 3 dechin dechin   4096 11月 18  2020 ./
drwxrwxr-x 3 dechin dechin   4096 7月  12 10:19 ../
-rw-r--r-- 1 dechin dechin 519979 11月 18  2020 ref.pdb
drwxr-xr-x 2 dechin dechin   4096 11月 18  2020 .solutions/
-rw-r--r-- 1 dechin dechin 117620 11月 18  2020 traj-broken.xtc
-rw-r--r-- 1 dechin dechin 116036 11月 18  2020 traj-whole.xtc

在具备了轨迹文件之后，就可以按照官方的文档示例来编写plumed.dat配置输入文件，用于定义需要计算的内容，比如这里定义的是计算1号原子和2号原子之间的距离，并每10个step将其写入到名为colvar的文件下：

(base) dechin@ubuntu2004:~/projects/gitlab/dechin/src/plumed/examples/trieste-1$ cat plumed.dat 
d1: DISTANCE ATOMS=1,2 
PRINT ARG=d1 FILE=colvar STRIDE=10

执行plumed指令得到的结果如下：

(base) dechin@ubuntu2004:~/projects/gitlab/dechin/src/plumed/examples/trieste-1$ plumed driver --mf_xtc traj-whole.xtc --plumed plumed.dat 

DRIVER: Found molfile format trajectory xtc with name traj-whole.xtc
PLUMED: PLUMED is starting
PLUMED: Version: 2.7.1 (git: Unknown) compiled on Jul 12 2021 at 09:24:30
PLUMED: Please cite these papers when using PLUMED [1][2]
PLUMED: For further information see the PLUMED web page at http://www.plumed.org
PLUMED: Root: /usr/local/lib/plumed
PLUMED: For installed feature, see /usr/local/lib/plumed/src/config/config.txt
PLUMED: Molecular dynamics engine: driver
PLUMED: Precision of reals: 8
PLUMED: Running over 1 node
PLUMED: Number of threads: 1
PLUMED: Cache line size: 512
PLUMED: Number of atoms: 6580
PLUMED: File suffix: 
PLUMED: FILE: plumed.dat
PLUMED: Action DISTANCE
PLUMED:   with label d1
PLUMED:   between atoms 1 2
PLUMED:   using periodic boundary conditions
PLUMED: Action PRINT
PLUMED:   with label @1
PLUMED:   with stride 10
PLUMED:   with arguments d1
PLUMED:   on file colvar
PLUMED:   with format  %f
PLUMED: END FILE: plumed.dat
PLUMED: Timestep: 1.000000
PLUMED: KbT has not been set by the MD engine
PLUMED: It should be set by hand where needed
PLUMED: Relevant bibliography:
PLUMED:   [1] The PLUMED consortium, Nat. Methods 16, 670 (2019)
PLUMED:   [2] Tribello, Bonomi, Branduardi, Camilloni, and Bussi, Comput. Phys. Commun. 185, 604 (2014)
PLUMED: Please read and cite where appropriate!
PLUMED: Finished setup
PLUMED:                                               Cycles        Total      Average      Minimum      Maximum
PLUMED:                                                    1     0.003111     0.003111     0.003111     0.003111
PLUMED: 1 Prepare dependencies                             5     0.000005     0.000001     0.000000     0.000003
PLUMED: 2 Sharing data                                     1     0.001410     0.001410     0.001410     0.001410
PLUMED: 3 Waiting for data                                 1     0.000022     0.000022     0.000022     0.000022
PLUMED: 4 Calculating (forward loop)                       1     0.000017     0.000017     0.000017     0.000017
PLUMED: 5 Applying (backward loop)                         1     0.000047     0.000047     0.000047     0.000047
PLUMED: 6 Update                                           1     0.000052     0.000052     0.000052     0.000052

在这个给定的轨迹文件中其实包含了5个step，但是由于我们设置了每10个step打印一次，因此最终存储到colvar中的数据只有1个值：

(base) dechin@ubuntu2004:~/projects/gitlab/dechin/src/plumed/examples/trieste-1$ ll
总用量 760
drwxr-xr-x 3 dechin dechin   4096 7月  12 10:22 ./
drwxrwxr-x 3 dechin dechin   4096 7月  12 10:19 ../
-rw-rw-r-- 1 dechin dechin     37 7月  12 10:22 colvar
-rw-rw-r-- 1 dechin dechin     60 7月  12 10:21 plumed.dat
-rw-r--r-- 1 dechin dechin 519979 11月 18  2020 ref.pdb
drwxr-xr-x 2 dechin dechin   4096 11月 18  2020 .solutions/
-rw-r--r-- 1 dechin dechin 117620 11月 18  2020 traj-broken.xtc
-rw-r--r-- 1 dechin dechin 116036 11月 18  2020 traj-whole.xtc
(base) dechin@ubuntu2004:~/projects/gitlab/dechin/src/plumed/examples/trieste-1$ cat colvar 
#! FIELDS time d1
 0.000000 0.095760

但是如果我们把上面的存储步长STRIDE修改为1的话，得到的colvar就是如下所示的结果：

(base) dechin@ubuntu2004:~/projects/gitlab/dechin/src/plumed/examples/trieste-1$ cat colvar 
#! FIELDS time d1
 0.000000 0.095760
 1.000000 0.096845
 2.000000 0.095885
 3.000000 0.096010
 4.000000 0.095321

这里是将5个step中的所有对应的原子距离都返回了出来。

轨迹文件的格式转换

上面用到的xtc文件是GROMACS生成的轨迹文件扩展名，由于是二进制文件并不方便读取，这里我们可以将其转换成hdf5的格式，然后就可以用python直接来读取其中的数据。这里我们需要借助于mdtraj这个工具。

用pip安装mdtraj

经过尝试用conda来装失败了，因此建议直接使用pip来进行安装：

(base) dechin@ubuntu2004:~/projects/gitlab/dechin/src/plumed/examples/trieste-1$ python3.9 -m pip install mdtraj
Collecting mdtraj
  Downloading mdtraj-1.9.6-cp39-cp39-manylinux_2_5_x86_64.manylinux1_x86_64.whl (6.0 MB)
     |████████████████████████████████| 6.0 MB 1.1 MB/s 
Collecting astunparse
  Using cached astunparse-1.6.3-py2.py3-none-any.whl (12 kB)
Requirement already satisfied: numpy>=1.6 in /home/dechin/anaconda3/envs/AmberTools21/lib/python3.9/site-packages (from mdtraj) (1.21.0)
Requirement already satisfied: scipy in /home/dechin/anaconda3/envs/AmberTools21/lib/python3.9/site-packages (from mdtraj) (1.7.0)
Requirement already satisfied: pyparsing in /home/dechin/anaconda3/envs/AmberTools21/lib/python3.9/site-packages (from mdtraj) (2.4.7)
Requirement already satisfied: six<2.0,>=1.6.1 in /home/dechin/anaconda3/envs/AmberTools21/lib/python3.9/site-packages (from astunparse->mdtraj) (1.16.0)
Requirement already satisfied: wheel<1.0,>=0.23.0 in /home/dechin/anaconda3/envs/AmberTools21/lib/python3.9/site-packages (from astunparse->mdtraj) (0.36.2)
Installing collected packages: astunparse, mdtraj
Successfully installed astunparse-1.6.3 mdtraj-1.9.6

安装完成后会生成一个名为mdconvert的可执行文件，一般情况是在系统路径下可以直接指令调用的，下述是mdconvert的帮助文档：

$ mdconvert -h
usage: mdconvert [-h] 
                 -o OUTPUT         # 必需，指定输出文件 
                 [-c CHUNK]        # 可选，指定一次读入内存的帧数，默认1000.
                 [-f]              # 可选，如果输出文件已存在，则强制覆盖 
                 [-s STRIDE]       # 可选，只加载 every stride-th frame 
                 [-i INDEX]        # 可选，更加灵活地指定加载特定的帧(Python格式)，eg： -i N; -i -1; -i N:M 
                 [-a ATOM_INDICES] # 可选，通过文件指定只加载特定的某些原子，文件内容为原子序号(从0开始计数)
                 [-t TOPOLOGY]     # 可选，指定拓扑文件(PDB/prmtop 格式)，一旦指定即可将轨迹输出为pdb格式。
                 input [input ...] # 必需，指定输入文件(可以是多条轨迹)

通过mdconvert，可以将上面的案例中所提到的xtc文件配合pdb文件转化成hdf5格式的文件：

(base) dechin@ubuntu2004:~/projects/gitlab/dechin/src/plumed/examples/trieste-1$ mdconvert -o test.h5 -t ref.pdb traj-whole.xtc 
Warning: 'step' data from input file(s) will be discarded. output format only supports fields: 'xyz', 'time', 'cell_lengths', 'cell_angles', 'velocities', 'kineticEnergy', 'potentialEnergy', 'temperature', 'lambda', 'topology'
converted 5 frames, 6580 atoms

而读取hdf5的文件，可以用python中的h5py来实现，没有安装h5py的可以通过pip简单的安装一下：

(base) dechin@ubuntu2004:~/projects/gitlab/dechin/src/plumed/examples/trieste-1$ python3 -m pip install h5py==2.9
Collecting h5py==2.9
  Downloading h5py-2.9.0-cp38-cp38-manylinux1_x86_64.whl (2.8 MB)
     |████████████████████████████████| 2.8 MB 997 kB/s 
Requirement already satisfied: numpy>=1.7 in /home/dechin/anaconda3/lib/python3.8/site-packages (from h5py==2.9) (1.20.2)
Requirement already satisfied: six in /home/dechin/.local/lib/python3.8/site-packages (from h5py==2.9) (1.16.0)
Installing collected packages: h5py
  Attempting uninstall: h5py
    Found existing installation: h5py 3.2.1
    Uninstalling h5py-3.2.1:
      Successfully uninstalled h5py-3.2.1
Successfully installed h5py-2.10.0

具体的h5py的使用方法这里就不进行展开了，可以参考官方的使用文档。

总结概要

本文作为一个入门级的文章，主要介绍了分子动力学模拟中增强采样的基本概念与相应工具的安装和使用。PLUMED是业界比较出名的一款增强采样开源软件，能够对接多个分子动力学模拟软件，如GROMACS等，并利用这些软件生成的路径信息来进行采样。并且为了可以在python上也能看到路径信息等重要数据，可以考虑使用mdconvert将路径数据转化成python上常用的hdf5格式并用h5py进行读写。