Self-learning Monte Carlo Method

YN, H. Shen, Y. Qi, J. Liu, and L. Fu“Self-learning Monte Carlo method: Continuous-time algorithm”,
Physical Review B 96, 161102(R) (2017) Editors’ Suggestion
原子力機構

システム計算科学センター

永井佑紀

Outline
イントロダクション
自己学習モンテカルロ法
いくつかの例：古典スピン系、二重交換模型
連続時間量子モンテカルロ法への適用
まとめ
YN, H. Shen, Y. Qi, J. Liu, and L. Fu“Self-learning Monte Carlo method: Continuous-time algorithm”,
Physical Review B 96, 161102(R) (2017) Editors’ Suggestion

マサチューセッツ工科大学(MIT)物理学科L. Fuグループとの共同研究

機械学習とは
猫画像その他関数f(x)y猫?

猫です

yes or no

学習プロセスx

大量の入力データxを用いて、f(x)=yを満たす関数f(x)を決める

例：囲碁棋譜データ→勝利の方程式
囲碁プロ棋士に勝利

自己学習モンテカルロ法

スピン配置を使って高速に

ボルツマン重みシミュレーション

分配関数と物理量物理量の期待値<A>を計算したい
平衡統計力学：分配関数が計算できればよい<A>=Tr [A exp(-βH)]/Z
Z = Tr [exp(-βH)]
経路積分表示なら経路積分表示でも何でも、分配関数の計算が重要

どのくらい大変かトレースTrの空間のサイズは？

量子多体問題の場合
固体物理学：原子を飛び移る電子がどう振る舞うか
銅酸化物高温超伝導体 etc.
有効模型：ハバード模型
upスピン いるorいない:2通り

x

downスピン いるorいない:2通り
4通り

1サイト

4 x 4 = 16通り
2サイト

4 x 4 x 4 = 64通り
3サイト

Nサイト=4N通り
18サイト

膨大な自由度

418=68719476736~687億次元

量子多体問題の場合
固体物理学：原子を飛び移る電子がどう振る舞うか条件を課す：粒子数保存、スピン保存→ブロック対角化できる
3サイト、upスピン1個、downスピン1個

3C1

x3C1

=9通り

Nサイト、upスピンn個、downスピンm個
xNCm
NCn
18サイト、upスピンdownスピンともに9個
18C9 x18C9=48620 x 48620 = 2363904400=23億次元

次元が大きすぎる！→モンテカルロ法

古典統計力学での例
分配関数

Aに関する期待値
可能な全ての配置を足す必要がある

古典スピン系（Ising模型の場合）それぞれの配置iは異なるエネルギーEiを持つ
E1

E2

E3
エネルギーが低い方が寄与が大きい→被積分関数が局在しているエネルギーが低い配置がよく出る偏ったサイコロでサンプリングすれば計算できる→重みつきモンテカルロ積分
偏ったサイコロはどうやって設計するか→マルコフ連鎖で可能自己学習モンテカルロ法

マルコフ連鎖モンテカルロ法
求める確率分布を均衡分布として持つマルコフ連鎖を生成→確率分布のサンプリングを行う

0.4

局在した関数π(x)を持つ高次元積分

π(x)

0.35

0.3

0.25

0.2

0.15

0.1

0.05

0-3-2-1x0

確率分布π(x)のマルコフ連鎖 →123

高精度に計算できる

非常に広い範囲で用いられる
ベイズ推定、量子化学、多体問題、宇宙・素粒子（格子QCD） etc.
例：古典統計力学
分配関数
物理量Aの統計平均

重みwi = exp(-βEi)→確率分布
高速化は様々な分野で恩恵が得られる

マルコフ連鎖モンテカルロ法
偏ったサイコロπ(x)を使えばただ足しあげるだけで済む

ランダムに配置xiを更新していく

w(x1)

w(x2)

w(x3)

x1=(1,-1,1,1)x2=(1,-1,-1,1)x3=(1,-1,-1,-1)

x1 x2 x3 x4 … xn

受け入れ確率の表式

1. 候補となる配置を生成する
どんな配置にするか？

2.その候補の受け入れ確率αを計算する 高い受け入れ確率が欲しい

メトロポリスヘイスティング法

従来のモンテカルロ法
配置の更新には2通りの方法がある

w(x1)

w(x2)

x1=(1,-1,1,1)x2=(1,-1,-1,1)

1. ローカルアップデート
局所的に一部だけ変化させる

長所：どのマルコフ連鎖モンテカルロ法でも可能
短所：変化しにくい（自己相関長が長い）ため、独立した配置が中々得られない

2. グローバルアップデート

高次元空間の探索が終わらない

w(x1)

w(x2)

x1=(1,-1,1,1)x2=(-1,1,-1,-1)

モデルに最適化された配置更新方法
長所：有効な方法が見つかればとても効率的

短所：モデルに強く依存する、探すのに職人的技巧が必要、存在するか非自明似たモデルでも使えない場合あり。見つけるだけで論文になる

自己学習モンテカルロ法
3. 自己学習モンテカルロ法(SLMC法)“モデルに最適化された配置更新方法”を学習で自分で見つける

確率分布：P(x)=w(x)0.4

Peﬀ(x)

0.35
0.4

0.3

P(x)

0.350.3Heﬀ(x)を見つけ出す

0.25

0.25

0.2

0.2

0.15
0.15

0.1
0.1

0.05

0-3-2-1x0xweff(x) = exp[-βHeff(x)]
0.05

0

1

2

3

例：w(x) = exp[-βH(x)]学習データ：配置{xi}
教師（正解）：重み{w(xi)}-3-2-10123配置の提案を有効模型Heﬀを使って行う
受け入れ確率α

完璧な模型なら100%自己学習モンテカルロ法の計算の流れ(a) トライアルシミュレーション(local update)(b)

本シミュレーション（SLMC法）

…

Heff

H

次配置提案採択する？

cj候補ci
機械学習

cj

Heff

1.トライアルシミュレーション

特徴：学習データを取りたいだけ取れる

1-1. 通常通り計算を行い、配置(学習データ)と重み(正解)を集める
1-2. 正解を再現するような有効模型Heﬀを、学習によって決定する

2.本シミュレーション
2-1. ある配置xiから、Heﬀを用いたマルコフ連鎖モンテカルロ法を行い、ほとんど相関のない配置cjを得る
2-2. Hを用いて受け入れ確率αを計算する

統計的に元シミュレーションに厳密に一致する

Heﬀ(x)を見つけ出す

SLMCは統計的に厳密に等しい
アクセプト率の計算はオリジナル模型を使用 -> 統計的に厳密な結果
1. 有効模型の導入

メトロポリスヘイスティング法によるアクセプト率の表式

Heff

2.有効模型の詳細釣り合い条件
ciからcjへの飛び移り確率を提案確率

3.とみなすciからcjへと飛び移る確率4.アクセプト率の表式

pprop: 提案確率
pacc: アクセプト確率

5. Heff を用いてアップデートを行い、相関のない配置を提案する

c1

c2

c3

c4

c5

c6

有効模型構築
トライアルシミュレーションで随時生成
サイズ4x1=(1,-1,1,1)x2=(1,-1,-1,1)

w(x1)

w(x2)

x3=(1,-1,-1,-1) ←input dataw(x3)←output data
サイズ1x1=(1,-1,1,1)f(x)f(x)=-βHeff(x)log w(x1)線形回帰ならHeff(x) = Σi ai gi(x)+bgi(x)の形を仮定しaiを決定する
学習データは好きなだけ増やせるいくつかの例：
古典スピン系、二重交換模型

古典スピン系の場合
二次元正方格子上の古典スピン系

J. Liu, Y. Qi, Z. Y. Meng, and L. Fu, Phys. Rev. B 95, 041101(R) (2017)

K=0では、この模型はIsing模型となる強磁性相互作用： J>0, K>0

J
K→Wolﬀ法という有効なグローバルアップデー
ト法がある

相転移が存在する
高温領域：常磁性相
低温領域：強磁性相

有効なグローバルアップデート法がない
臨界点近傍で自己相関時間が非常に長い

SLMC法を使えば、臨界点近傍での自己相関時間を劇的に短縮することができる

古典スピン系の場合
二次元正方格子上の古典スピン系 J. Liu, Y. Qi, Z. Y. Meng, and L. Fu, Phys. Rev. B 95, 041101(R) (2017)

オリジナル模型

有効模型

二体相互作用のみ：
Wolﬀグローバルアップデート法が使える４体相互作用:
有効なグローバルアップデートがない

~
J1

J
K

配置と重みのデータを集める~
J2線形回帰Log w=-βH(Si)Log weff=-βHeff(Si)古典スピン系の場合
J. Liu, Y. Qi, Z. Y. Meng, and L. Fu, Phys. Rev. B 95, 041101(R) (2017)

オリジナル模型

有効模型

自己相関関数の減少の様子

C1:最近接スピンスピン相関関数

J1 だけを用いて元の模型の重みを再現できる

量子系
フェルミオン系

J. Liu, H. Shen, Y. Qi, Z. Y. Meng, and L. Fu, Phys. Rev. B 95, 241104(R)(2017)

重みとして行列の行列式を計算しなければならない
φ:配置重い三次元立方格子上での二重交換模型
遍歴電子がスピンを介してRKKY相互作用を発現させる
ci: サイトiでのフェルミオン
Si: サイトiでの古典スピン

J/t=16 かつ1/4フィリングで、この模型は常磁性・強磁性転移を起こす（Tc〜0.12t）Tc近傍でのシミュレーションがしたい
有効模型：二体の局所スピン間の相互作用の模型
注：フェルミオン自由度を持たない

重みを再現するHeffを学習する

量子系
J. Liu, H. Shen, Y. Qi, Z. Y. Meng, and L. Fu, Phys. Rev. B 95, 241104(R)(2017)

二重交換模型：フェルミオン＋古典スピン有効模型

行列式

ただの数振動関数-> RKKY 相互作用

二重交換模型のフェルミオン自由度を消去

自己相関時間ほぼ同じ分布関数

L=8 (8x8x8 cubic)
103 倍の高速化!

->RKKY相互作用が出るRKKYを自動的に”導出”システムサイズ連続時間量子モンテカルロ法への適用

連続時間量子モンテカルロ法
連続時間量子モンテカルロ法（CTQMC）は動的平均場理論（DMFT）における不純物ソルバーとして使われることが多い
鈴木トロッター分解がいらない

数値的に厳密な計算ができる

DMFT (and its extensions): 銅酸化物超伝導体、鉄系超伝導体etc
DMFT+density functional theory (DFT):強相関系の第一原理計算
不純物アンダーソン模型

自由な電子と、局在した電子

連続時間量子モンテカルロ法
分配関数を H1で摂動展開する

ふた通りの展開方法: U or V as H1

E. Gull et al., Rev. Mod. Phys. 83, 349 (2011)

continuous-time auxiliary-field quantum Monte Carlo method (CTAUX)

配置0
イジング的な擬スピンを導入する

β τ連続時間量子モンテカルロ法
不純物アンダーソン模型

continuous-time auxiliary-field method (CTAUX)

イジング的な擬スピン分配関数

配置0ファインマンダイアグラムを重みとみなしてモンテカルロ積分β τβ τ0

連続虚時間上のスピンの数はシミュレーション中に変化する

有効模型は何か？

連続時間量子モンテカルロ法
YN, H. Shen, Y. Qi, J. Liu and L. Fu, Phys. Rev. B 96, 161102(R) (2017)

元の模型

有効模型

虚時間上で長距離の二体スピンスピン相互作用異なる展開次数（ファインマンダイヤグ
重みの計算負荷が大きい
トレーニングデータ
インプット: 連続虚時間上のスピンの配置スピン数はシミュレーション中に変化ラム）の重みを30個程度のパラメータで再現ができる！
その他面白い性質
1. 最低次のスピンスピン相互作用のみ

アウトプット: 行列式を計算し2. 異なる展開次数の重みが同一の関数で表現できるて得られる重み3. 1/n因子は原子極限から導出されており、これが無いと模型の精度が大きく落ちる

連続時間量子モンテカルロ法(CTAUX)CTAUXの重みの計算

有効模型：虚時間上で二体相互作用するスピン

YN, H. Shen, Y. Qi, J. Liu and L. Fu, Phys. Rev. B 96, 161102(R) (2017)

連続変数→多項式展開で有限次数で打ち切り行列式の計算が必要

全部で30個ほどのパラメータを線形回帰で決定重みの分布関数

有効模型の長距離相互作用
5(J(τ))/|U| with U=2(J(τ))/|U| with U=5(L(τ)+b0)/|U| with U=2(L(τ)+b0)/|U| with U=54
32
10-1-2-3-400.20.4

τ/β0.60.81

有効模型は、出てくるファインマンダイヤグラムの全ての重要な寄与を記述している

計算コストOriginal without fast updates:
O(<n>3 τori)
Original with fast updates:
O(<n>2 τori)
Self-learning MC:
O((<n>3+mc <n>τori)τSL)

τori ∝ U4β3

Markov Chain Monte Carlo steps
Original CTAUX

Self-learning CTAUX:Local fast update O(n2):Local determinant update O(n3):Local self-learning update O(n)自己学習モンテカルロ法は自己相関時間を減らす

まとめ自己学習モンテカルロ法（SLMC）“モデルに最適化された配置更新方法”を学習で自分で見つける
YN, H. Shen, Y. Qi, J. Liu and L. Fu, Phys. Rev. B 96, 161102(R) (2017)
Trial simulations

Main simulations（SLMC）…
Heff

Machine learning

xi

Heff

H

ProposeAccept？

proposal xjxj

ポイント：次配置提案に機械学習が使われる->物理量計算はオリジナルを使う->オリジナル手法と統計的に厳密に等しい
SLMCが有効な時
1. 自己相関時間が長い (e.g. 臨界点近傍)
2. 似た系が有効なグローバルアップデートを持っている
3. 重みの計算が大変（フェルミオン系等の量子系）