誤差の評価
ここでは誤差について簡単にまとめています。特にモンテカルロ計算での統計誤差の付け方を説明しています。
素粒子やハドロン物理分野で統計誤差を付けるとき、その多くの場合でデータの分布は正規分布に従います。
そのため、誤差は正規分布に従って 1σ信頼区間や2σ信頼区間を付けることになります。
もちろん、非常にまれですが正規分布にならない場合もあり得るので注意は必要です(特にコーシー分布)。
他分野では P値を計算することもありますが、本研究分野では基本的にありません。
とは言え無いとは限らないので、その際は多重性に注意して解析してください(そうしないと意図していないのに P値ハッキングをしてしまいます)。
例えば、AMED(国立研究開発法人日本医療研究開発機構)が公開している資料「多重性の問題:研究計画の重要性」を確認すると良いでしょう。
ジャックナイフ法
モンテカルロ計算はマルコフ鎖上で配位を生成し、それを用いて物理量の期待値を計算します。
しかしその性質上、ある配位はその前後の配位と相関を持つため、単純に統計誤差を計算すると誤差を過小評価します。
これは、マルコフ鎖上である程度離して(ある程度の配位を捨てながら一定間隔で配位を保存して)期待値を計算することで抑えられます。
しかし、計算コストの都合上、多くの場合で相関が残り誤差を過小評価します。
この相関(自己相関と呼びますが)を見積もる方法の一つにジャックナイフ法があります。
以下では、その計算過程を見ていきましょう。
まず、配位が N 個(十分な個数)あるとし、これをいくつかの集合に分割します。
この分割された各集合をビンと呼び、ビンの要素(つまりそこに含まれている配位の数)はすべてのビンで同じとします。
ビンの数は M とすると、ビンの中の配位の数は L = N/M です。
次に、ビンを一つ取り除いた M-1 個のビンを用いて統計平均を計算します。
これを <O>_b = \sum_k O_k/(N-L) と書きましょう。
ここで、\sum は N-L 個の配位を用いて和を取ることを意味し(一つビンを取り除いているので)、O_k は計算したい物理量の各配位での値とします。
この平均値を用いて、
<O>=\sum_b O_b /M, δ = \sqrt[(L-1)(<O>^2-<O^2>)]
と期待値 <O> とその統計誤差 δ を計算します。
ここで \sum_b は、一つビンを取り除いた集合の全ての組み合わせについての和の意味です。
後はビンの数を変えてその誤差の振る舞いを確認して、誤差の増減がほぼ無く一定になった値を最終的な統計誤差とします。
適当に書いたコードも載せておきますが、バグチェックはしていないので参考までに。
使い方はコードを読んで確認してください。
この方法の詳細は、
[1] 格子上の場の理論、青木真也 著、丸善出版
[2] ゼロからできるMCMC マルコフ連鎖モンテカルロ法の実践的入門、花田政範・松浦壮 著、講談社
等を参照すると良いでしょう。
他の手法
ジャックナイフ法以外でモンテカルロ法の結果に統計誤差を付けるとしたらどうすればよいでしょうか?
一つの方法は独立なマルコフ鎖を複数用意してそれぞれから一つづつ配位を取る方法です。
この場合、統計誤差は通常の手法で付けることができます。
この手法の問題はいくつも並列して計算しないと時間がかかること、熱化などで捨てる配位まで考えるとコストがかかる点でしょうか。
時間とコストが問題にならない場合はこの手法でもよいでしょう。
この場合は基本的に(十分初期状態から話してあれば)自己相関が無いので誤差がつけやすいですね。