副本交换动力学

许多分子模拟方案都需要对具有多个极小值的能量空间进行遍历采样，并且极小值之间的势垒在环境温度范围内难以跨越可访问的模拟时间范围。

这意味着结果与初始条件的选择相混淆，因为初始条件决定了模拟所探索的空间区域。

副本交换模拟通过在不同的系综中运行多个独立副本，并在这些系综之间定期交换副本的坐标来寻求增强这种情况下的采样。

典型地，我们构造一组副本，使得其中的一个极端是希望对其进行采样的系综，而另一极端是容易越过势垒的系综。

因此，如果尝试在系综之间交换副本会产生有效的系综

交换尝试成功的可能性足够高

副本在整个系综上的流动会产生良好的混合效果

那么所得的采样在统计上将是正确的，并且更接近遍历理想。

这怎么工作？它依赖于在系综之间进行副本的蒙特卡罗移动。在特定状态下观察到副本的概率取决于势能和温度。

观察到的概率范围服从正态分布（由于中心极限定理）。

如果我们选择两个系综，使它们的状态分布具有明显的重叠，那么当我们在一个系综中观察到一个状态时，很有可能在另一个系综中也观察到同样的状态。

根据同时在两个系综中观测到的概率来尝试交换在两个副本中找到的坐标，仅当随机数具有合适的值时，才接受交换。

如果正常完成交换，则可以实现细微的平衡；无论是否进行交换，两个系综中的采样都是正确的。

应该选择什么范围的系综？ REMD方法（也称为平行淬火）在一定温度范围内进行相同的模拟。选择最高温度，以便越过势垒。

接受交换后，副本会在温度空间中上下徘徊；在所有温度下都可能会越过势垒，但越高的温度下会越容易。如果这些状态“属于”相应的系综，则会下降到较低的温度。这样，可以提高采样的遍历性。

副本交换分子动力学(REMD)是一种增强采样的技术，采用的方法是在不同温度下对具有相似势能的体系进行采样。通过这种方法，体系可能会越过势垒，从而探索新的构象空间。确定了多肽及其周围环境后，就需要确定待采样温度空间的范围、使用的处理器数目，以及模拟的时间长度。体系、副本数目、温度空间的范围和温度分布决定了副本之间的平均交换概率。这些值对所有副本都应该是相同的。为此，如果假定自由能空间中没有已知的瓶颈部分，那么预期体系的势能会随温度的增加而增加，因此副本的温度分布应服从指数分布。

有充分的证据表明，对多肽/蛋白质体系，副本尝试交换的时间间隔不应小于1 ps。这可以决定在两次交换之间，一个副本探索相应构象空间所用的平均时间，这个时间比实际的交换概率还重要。在GROMACS的实现中，特定一对副本的交换尝试会隔次进行，因为奇数对和偶数对尝试交换会交替进行。注意，模拟是否“成功”完全取决于要解决的问题和判断成功的标准。虽然REMD模拟可以增加采样量，但并不能提供最终答案。这一点应该牢记在心。

定义体系，例如：多肽+溶剂(隐式或显式)

根据可用的处理器数目和要采样的温度范围(实际上它们的相关性非常强)，来选择温度分布。其中 k和T可以微调，以获得合理的温度间隔，保证能够进行交换。指数形式可以保证温度间隔随温度升高而增大。这是必要的，因为总能量的分布随着温度的增加而增加，因此交换概率也随之增加。保持交换概率不随温度变化。

获得温度分布后，在N个温度下分别对体系进行平衡(每个温度使用单独的.mdp文件)。然后根据平衡结构创建一系列(N个)运行输入文件(.tpr)，仍要使用不同的.mdp文件生成不同的.tpr文件。

然后，运行短时间的REMD模拟以获得交换概率的估计值(大约100 ps就可以获得良好的估计值)，如果得到的结果与所需值偏离较大，试着修改温度。典型情况下，交换概率取0.2到0.3是合适的。实际上，比交换概率本身的值更重要的是副本在给定温度下的停留时间。它的值为交换概率与交换频率的乘积。例如，交换概率为0.2，每隔2 ps进行一次交换尝试，可以得出相应温度下的平均时间为10 ps。另外还需要考虑的是交换方式。每次尝试时，是尝试交换所有的副本对，还是随机选择一对副本进行交换。对后一种情况，还应考虑一点，必须将副本的停留时间乘以N-1，即交换对的数目。

每个温度下的初始构象可以相同，也可以不同。如何选择取决于运行REMD的原因。

网上有生成T-REMD温度的源码，它是基于溶剂化蛋白已知的能量分布计算的。输入蛋白质的原子数，水的原子数，温度范围和所需的交换概率，就会生成所需的温度。

副本交换的运行步数，即要达到所需的交换尝试周期所需的步数(由mdrun选项-replex指定)。

mpirun -np x gmx_mpi mdrun -s topol -replex m -multidir <目录名称>

注意：

如果多重模拟中的模拟顺序很重要，那么当将其名称提供给-multidir时，需要对它们的名称进行排序

应谨慎使用支持通配符的shell，因为可能出现这种情况：dir1 dir10 dir11 ... dir2...

组截止模式的gmx mdrun -table仅在此模式下有效

mpirun -np 32 gmx_mpi mdrun -replex 1000 -multidir a b c d -v

使用4个模拟在32个进程上进行一个多重模拟。输入和输出文件位于目录a，b，c和d中。

mpirun -np 32 gmx_mpi mdrun -replex 1000 -multidir a b c d -gputasks 0000000011111111 -v

进行与之前相同的多重模拟。在具有两个物理节点、每个节点两个GPU的机器上，每个节点有16个MPI rank，而每个模拟用8个MPI rank。

在节点上执行PP工作的16个MPI rank被映射到ID为0和1的GPU，即使它们来自多个模拟。

它们按照指示的顺序进行映射，以便每个模拟的PP rank使用单个GPU。但是，使用0101010101010101的顺序，运行速度可能更快。

副本交换分子动力学（Replica Exchange MolecularDynamics，REMD）是一种加速采样的方法，可用于任何类型的模拟，特别适用于不同构象之间存在较高能量势垒的情况。它对同一系统处于不同温度下的多个副本进行模拟，并按一定的时间间隔随机地交换两个副本的完整状态，交换概率为：

其中T1 和T2为参考温度，U1和U2分别为副本 1 和 2 的瞬时势能。交换后，使用（T1/T2)±0.5 对速度进行缩放，并在下一步进行邻区搜索。这种方法结合了最高温度的快速采样和频繁跨垒，以及所有不同温度下正确的 Boltzmann 采样。因为概率随着温差的增加而迅速下降，我们只会尝试与相邻的温度进行交换。在一步中不应该尝试对所有可能的配对副本进行交换。例如，如果交换副本 1 和2，交换副本 2 和 3 的机会不仅取决于副本2和3的能量，还取决于副本 1 的能量。在 GROMACS 中，对这个问题的解决方法是，在 奇数次尝试时试着交换所有的 奇数副本对，在 偶数次尝试时试着交换所有的偶数副本对。如果我们有四个按温度排序的副本：0，1，2 和 3，并且每 1000 步尝试交换一次，那么将会在 1000，3000 等步尝试交换 0-1 和 2-3，在 2000，4000 等步尝试交换 1-2。

应该如何选择温度呢? 能量差可以写成:

其中Ndf为一个副本的总自由度数，c对简谐势为 1，对蛋白质/水系统为 2 左右。如果 T2 = (1 + ε)T1则概率变为：

因此对于大小为e−2 ≈0.135 的概率，可得到 𝜖≈2/√cNdf。所有键都受到约束时，Ndf≈2 Natoms，因此，对于c=2，应选择 为 1/√ Natoms。但是，当使用压力耦合时存在一个问题。在较高的温度下，密度会降低，从而导致能量增高，这些我们也应该考虑到。GROMACS 网站提供了一个所谓的 REMD计算器，可以根据你输入的温度范围和原子数给出一组建议温度。

Okabe等人提出了等压等温系综的REMD扩展. 他们将交换概率修改为:

其中P1和P2为各自的参考压力，V1和V2为模拟中各自的瞬时体积。在大多数情况下，体积的差别非常小，第二项可以忽略不计。只有当P1和P2之间相差很大或存在相变时，第二项才起作用。

GROMACS 还支持哈密顿量副本交换。在哈密顿量副本交换中，每个副本的哈密顿量都不同，是由模拟指定的自由能路径定义的。用来维持正确系综概率的交换概率为

也可以同时进行哈密顿量副本交换和温度副本交换, 使用的接受标准为:

mdrun程序的包含了所有不同副本交换方法的选项. 由于算法中内在的并行化, 只有安装了MPI时候才能运行副本交换模拟. 为提高效率, 每个副本可以在一个独立的队列上运行. 想了解如何使用这些多节点功能, 可以参考mdrun的帮助文档.