加权直方图分析方法 (WHAM)

这是一个用 Rust 编写的快速加权直方图分析方法实现。它允许从雨伞采样模拟中计算多维自由能轮廓。有关方法的更多详细信息，建议参考 Roux, B. (1995). 使用计算机模拟计算平均力势，CPC，91(1)，275-282。

功能

• 快速，特别是对于小系统
• 多线程（自动在所有可用核心上运行）
• 多维（任何数量的集体变量都是可能的）
• 自相关分析以去除相关样本
• 通过重抽样进行错误分析
• 单元测试

安装

通过 cargo 从源代码安装

# cargo installation
curl -sSf https://static.rust-lang.org/rustup.sh | sh

cargo install wham

用法

wham 有一个方便的命令行界面。您可以使用 wham -h 查看所有选项

wham 1.1.0
D. Bauer <bauer@cbs.tu-darmstadt.de>
wham is a fast implementation of the weighted histogram analysis method (WHAM) written in Rust. It currently supports
potential of mean force (PMF) calculations in multiple dimensions at constant temperature.

Metadata file format:
    /path/to/timeseries_file1  x_1  x_2  x_N  fc_1  fc_2  fc_N
    /path/to/timeseries_file2  x_1  x_2  x_N  fc_1  fc_2  fc_N
    /path/to/timeseries_file3  x_1  x_2  x_N  fc_1  fc_2  fc_N
The first column is a path to a timeseries file _relative_ to the metadata file (see below). This is followed by the
position of the umbrella potential x in N dimensions and the force constant fc in each dimension. Lines starting with a
# are treated as comments and will not be parsed.

Timeseries file format:
    time  x_1  x_2  x_N
    time  x_1  x_2  x_N
    time  x_1  x_2  x_N
The first column will be ignored and is followed by N reaction coordinates x.

Shipped under the GPLv3 license.

USAGE:
    wham [FLAGS] [OPTIONS] --bins <BINS> --max <HIST_MAX> --file <METADATA> --min <HIST_MIN> --temperature <temperature>

FLAGS:
    -c, --cyclic     For periodic reaction coordinates. If this is set, the first and last coordinate bin in each
                     dimension are treated as neighbors for the bias calculation.
    -h, --help       Prints help information
    -g, --uncorr     Estimates statistical inefficiency of each timeseries via autocorrelation and removes correlated
                     samples (default is off).
    -V, --version    Prints version information
    -v, --verbose    Enables verbose output.

OPTIONS:
    -b, --bins <BINS>                  Number of histogram bins (comma separated).
        --bt <bootstrap>               Number of bayesian bootstrapping runs for error analysis by assigning random
                                       weights (defaults to 0).
        --seed <bootstrap_seed>        Random seed for bootstrapping runs.
        --convdt <convdt>              Performs WHAM for slices with the given delta in time and returns an output file
                                       for each slice. THis is useful to check the result for convergence. Example: with
                                       --convdt 100 and a timeseries ranging from 0-300, free energy surfaces for slices
                                       0-100, 0-200 and 0-300 will be given returned.
        --end <end>                    Skip rows in timeseries with an index larger than this value (defaults to 1e+20)
    -i, --iterations <ITERATIONS>      Stop WHAM after this many iterations without convergence (defaults to 100,000).
        --max <HIST_MAX>               Histogram maxima (comma separated). Also accepts "pi".
    -f, --file <METADATA>              Path to the metadata file.
        --min <HIST_MIN>               Histogram minima (comma separated for multiple dimensions). Also accepts "pi".
    -o, --output <output>              Free energy output file (defaults to wham.out).
        --start <start>                Skip rows in timeseries with an index smaller than this value (defaults to 0)
    -T, --temperature <temperature>    WHAM temperature in Kelvin.
    -t, --tolerance <TOLERANCE>        Abortion criteria for WHAM calculation. WHAM stops if abs(F_new - F_old) <
                                       tolerance (defaults to 0.000001).

示例

示例文件夹包含两个简单测试系统的输入和输出文件

• 1d_cyclic：真空中二氨酸的 Phi 拐角
• 2d_cyclic：同一系统的 Phi 和 psi 拐角

以下命令将运行二维示例（二氨酸的 phi 和 psi 角的模拟）并基于两个在 -3.14 到 3.14 范围内的集体变量计算自由能，每个维度有 100 个 bin 和周期性集体变量

wham --max 3.14,3.14 --min -3.14,-3.14 -T 300 --bins 100,100 --cyclic -f example/2d/metadata.dat       
> Supplied WHAM options: Metadata=example/2d/metadata.dat, hist_min=[-3.14, -3.14], hist_max=[3.14, 3.14], bins=[100, 100] verbose=false, tolerance=0.000001, iterations=100000, temperature=300, cyclic=true
> Reading input files.
> 625 windows, 624262 datapoints
> Iteration 10: dF=0.389367172324539
> Iteration 20: dF=0.21450559607810152
(...)
> Iteration 620: dF=0.0000005800554892309461
> Iteration 630: dF=0.00000047424278621817084
> Finished. Dumping final PMF
(... pmf dump ...)

收敛后，最终的偏差偏移（F）和自由能将被输出到 stdout，并写入输出文件。

输出文件包含每个 bin 的自由能和概率。概率被归一化，总和为 P=1.0，最小的自由能被设置为 0（其他自由能基于此）。

#coord1    coord2    Free Energy    +/-    Probability    +/-
-3.108600    	-3.108600    10.331716    0.000000    0.000095    0.000000
-3.045800    	-3.108600    8.893231    0.000000    0.000170    0.000000
-2.983000    	-3.108600    7.372765    0.000000    0.000312    0.000000
-2.920200    	-3.108600    6.207354    0.000000    0.000498    0.000000
-2.857400    	-3.108600    4.915298    0.000000    0.000836    0.000000
-2.794600    	-3.108600    3.644738    0.000000    0.001392    0.000000
-2.731800    	-3.108600    3.021743    0.000000    0.001787    0.000000
-2.669000    	-3.108600    2.827463    0.000000    0.001932    0.000000
-2.606200    	-3.108600    2.647531    0.000000    0.002076    0.000000
(...)

错误分析

WHAM可以使用贝叶斯自助法进行误差分析。每个模拟窗口被认为是一组独立的数据点。通过为每个窗口随机分配权重并计算概率N次，可以估计误差为N次自助运行之间的标准偏差。更多详情请参阅Van der Spoel, D. 等. (2010). g_wham—一种包括鲁棒误差和自相关估计的免费加权直方图分析实现，JCTC，6(12)，3713-3720。

在WHAM中执行贝叶斯自助法，使用-bt <RUNS>标志进行单个自助运行。箱概率和自由能的误差估计将以标准误差（SE）的形式在输出文件中的单独列（+/-）中给出。如果没有执行误差分析，这些列将被设置为0.0。

自相关分析

使用--uncorr标志，WHAM计算所有时间序列和所有集体变量的自相关时间tau。然后根据它们最高的自相关时间对这些时间序列进行过滤，以从数据集中去除相关样本。这减少了数据点的数量，但可以改善结果精度。

对于过滤，计算统计效率g： g = 1 + 2*tau，并且只使用时间序列的每g个元素进行去偏。该方法的更详细描述可参见Chodera, J.D. 等. (2007). 使用加权直方图分析法分析模拟和并行温度模拟，JCTC 3(1):26-41

许可证 & 引用

WHAM在GPL-3.0许可证下发布。请阅读此存储库中的LICENSE文件以获取更多信息。

此WHAM实现没有出版物。然而，有一个可引用的DOI。如果您使用此软件进行工作，请考虑引用它：Bauer, D.，WHAM - 使用Rust编写的有效加权直方图分析实现，Zenodo. https://doi.org/10.5281/zenodo.1488597

本工作中的一些部分，尤其是某些性能优化和I/O格式，受到了A. Grossfield的实现（Grossfield, A，WHAM：加权直方图分析法，http://membrane.urmc.rochester.edu/content/wham）的启发。