交換モンテカルロ法 `exchange`¶

exchange はレプリカ交換モンテカルロ法を用いてパラメータ探索を行う Algorithm です。

前準備¶

あらかじめ mpi4py をインストールしておく必要があります。:

python3 -m pip install mpi4py

入力パラメータ¶

サブセクション param と exchange を持ちます。

[`param`] セクション¶

initial_list

形式: 実数型のリスト。長さはdimensionの値と一致させます。

説明: パラメータの初期値。定義しなかった場合は一様ランダムに初期化されます。
unit_list

形式: 実数型のリスト。長さはdimensionの値と一致させます。

説明: 各パラメータの単位。
探索アルゴリズム中では、各パラメータをそれぞれこれらの値で割ることで、簡易的な無次元化・正規化を行います。定義しなかった場合にはすべての次元で 1.0 となります。
min_list

形式: 実数型のリスト。長さはdimensionの値と一致させます。

説明: パラメータが取りうる最小値。
max_list

形式: 実数型のリスト。長さはdimensionの値と一致させます。

説明: パラメータが取りうる最大値。

[`exchange`] セクション¶

numsteps

形式: 整数値。

説明: モンテカルロ更新を行う回数。
numsteps_exchange

形式: 整数値。

説明: 「温度」のレプリカ交換を行う頻度。この回数だけモンテカルロ更新を行ったらレプリカ交換を実行します。
Tmin

形式: 実数値。 (default: 0.1)

説明: 「温度」の最小値。
Tmax

形式: 実数値。 (default: 10.0)

説明: 「温度」の最大値。
Tlogspace

形式: 真偽値。 (default: true)

説明: 「温度」を各レプリカに割り当てる際に、対数空間で等分割するか否かを指定します。
true のときは対数空間で等分割します。

出力ファイル¶

`RANK/trial.txt`¶

各レプリカについて、モンテカルロサンプリングで提案されたパラメータと、対応する目的関数の値です。 1列目にステップ数、2列目にレプリカの温度、3列目に目的変数の値、4列目以降にパラメータが記載されます。

# step T fx z1 z2
0.004999999999999999 0.07830821484593968 3.682008067401509 3.9502750191292586
0.004999999999999999 0.0758494287185766 2.811346329442423 3.691101784194861
0.004999999999999999 0.08566823949124412 3.606664760390988 3.2093903670436497
0.004999999999999999 0.06273922648753057 4.330900869594549 4.311333132184154

`RANK/result.txt`¶

各レプリカについて、モンテカルロサンプリングで生成されたパラメータと、対応する目的関数の値です。 trial.txt と同一の書式です。

# step T fx z1 z2
0.004999999999999999 0.07830821484593968 3.682008067401509 3.9502750191292586
0.004999999999999999 0.07830821484593968 3.682008067401509 3.9502750191292586
0.004999999999999999 0.07830821484593968 3.682008067401509 3.9502750191292586
0.004999999999999999 0.06273922648753057 4.330900869594549 4.311333132184154

`best_result.txt`¶

サンプリングされた全データのうち、目的関数の値が最小となったパラメータと、対応する目的関数の値です。

nprocs = 4
rank = 2
step = 65
fx = 0.008233957976993406
z1 = 4.221129370933539
z2 = 5.139591716517661

アルゴリズム解説¶

マルコフ連鎖モンテカルロ法¶

モンテカルロ法（モンテカルロサンプリング）では、パラメータ空間中を動き回る walker \(\vec{x}\) を重み \(W(\vec{x})\) に従って確率的に動かすことで目的関数の最適化を行います。重み \(W(\vec{x})\) として、「温度」 \(T > 0\) を導入して \(W(\vec{x}) = e^{-f(\vec{x})/T}\) とすることが一般的です（ボルツマン重み）。ほとんどの場合において、 \(W\) に基づいて直接サンプリングする (walker を生成する) のは不可能なので、 walker を確率的に少しずつ動かすことで、頻度分布が \(W\) に従うように時系列 \(\{\vec{x}_t\}\) を生成します (マルコフ連鎖モンテカルロ法, MCMC)。 \(\vec{x}\) から \(\vec{x}'\) へ遷移する確率を \(p(\vec{x}' | \vec{x})\) とすると、

\[W(\vec{x}') = \sum_{\vec{x}} p(\vec{x}' | \vec{x}) W(\vec{x})\]

となるように \(p\) を定めれば時系列 \(\{\vec{x}_t\}\) の頻度分布が \(W(\vec{x})\) に収束することが示されます（釣り合い条件） 1 。実際の計算では、より強い制約である詳細釣り合い条件

\[p(\vec{x} | \vec{x}') W(\vec{x}') = W(\vec{x})p(\vec{x}' | \vec{x})\]

を課すことがほとんどです。両辺で \(vec{x}\) についての和を取ると釣り合い条件に帰着します。

\(p\) を求めるアルゴリズムはいくつか提案されていますが、 2DMAT では Metropolis-Hasting 法 (MH法) を用います。 MH 法では、遷移プロセスを提案プロセスと採択プロセスとに分割します。

提案確率 \(P(\vec{x} | \vec{x}_t)\) で候補点 \(\vec{x}\) を生成します
- 提案確率 \(P\) としては \(\vec{x}_t\) を中心とした一様分布やガウス関数などの扱いやすいものを利用します
提案された候補点 \(\vec{x}\) を採択確率 \(Q(\vec{x}, | \vec{x}_t)\) で受け入れ、 \(\vec{x}_{t+1} = \vec{x}\) とします - 受け入れなかった場合は \(\vec{x}_{t+1} = \vec{x}_t\) とします

採択確率 \(Q(\vec{x} | \vec{x}_t)\) は

\[Q(\vec{x} | \vec{x}_t) = \min\left[1, \frac{W(\vec{x})P(\vec{x}_t | \vec{x}) }{W(\vec{x}_t) P(\vec{x} | \vec{x}_t)} \right]\]

とします。この定義が詳細釣り合い条件を満たすことは、詳細釣り合いの式に代入することで簡単に確かめられます。特に、重みとしてボルツマン因子を、提案確率として対称なもの \(P(\vec{x} | \vec{x}_t) = P(\vec{x}_t | \vec{x})\) を用いたときには、

\[Q(\vec{x} | \vec{x}_t) = \min\left[1, \frac{W(\vec{x})}{W(\vec{x}_t)} \right] = \min\left[1, \exp\left(-\frac{f(\vec{x}) - f(\vec{x}_t)}{T}\right) \right]\]

という更に簡単な形になります。

\(\Delta f = f(\vec{x}) - f(\vec{x}_t)\) とおいて、 \(\Delta f \le 0\) のときに \(Q = 1\) となることを踏まえると、 MH 法によるMCMC は次のようになります。

現在地点の近くからランダムに次の座標の候補を選び、目的関数 \(f\) の値を調べる
\(\Delta f \le 0\) ならば（山を下る方向ならば）移動する
\(\Delta f > 0\) ならば採択確率 \(Q = e^{-\Delta f / T}\) で移動する
1-3 を適当な回数繰り返す

得られた時系列のうち、目的関数の値が一番小さいものを最適解とします。 3 番のプロセスのおかげで、 \(\Delta f \sim T\) ぐらいの山を乗り越えられるので、局所最適解にトラップされた場合にも脱出可能です。

レプリカ交換モンテカルロ法¶

モンテカルロ法による最適化では、温度 \(T\) は非常に重要なハイパーパラメータとなっています。モンテカルロ法では、温度 \(T\) 程度の山を乗り越えられますが、逆にそれ以上の深さの谷からは容易に脱出できません。そのため、局所解へのトラップを防ぐためには温度を上げる必要があります。一方で、 \(T\) よりも小さい谷は谷として見えなくなるため、得られる \(\min f(\vec{x})\) の精度も \(T\) 程度になり、精度を上げるためには温度を下げる必要があります。ここから、最適解を探すためには温度 \(T\) を注意深く決める必要があることがわかります。

この問題を解決する方法として、温度 \(T\) を固定せずに更新していくというものがあります。たとえば、焼きなまし法 (simulated annealing) では、温度をステップごとに徐々に下げていきます。焼戻し法 (simulated tempering) は、温度をハイパーパラメータではなく、サンプリングすべきパラメータとして扱い、（詳細）釣り合い条件を満たすように更新することで、加熱と冷却を実現します。温度を下げることで谷の詳細を調べ、温度を上げることで谷から脱出します。レプリカ交換モンテカルロ法 (replica exchange Monte Carlo) は焼戻し法を更に発展させた手法で、並列焼戻し法 (parallel tempering) とも呼ばれます。レプリカ交換モンテカルロ法では、レプリカと呼ばれる複数の系を、それぞれ異なる温度で並列にモンテカルロシミュレーションします。そして、ある一定間隔で、（詳細）釣り合い条件を満たすように他のレプリカと温度を交換します。焼戻し法と同様に、温度を上下することで谷を調べたり脱出したりするのですが、各温度点について、かならずレプリカのどれかが対応しているため、全体として特定の温度に偏ることがなくなります。また、一つのMPI プロセスに一つのレプリカを担当させることで簡単に並列化可能です。数多くのレプリカを用意することで温度間隔が狭まると、温度交換の採択率も上がるため、大規模並列計算に特に向いたアルゴリズムです。有限温度由来の「ぼやけ」がどうしても生まれるので、モンテカルロ法の結果を初期値として minsearch をするのがおすすめです。

脚注

1: 正確には、収束のためには非周期性とエルゴード性も必要です。

2DMAT 1.0.1 ドキュメント