57
ＭＣＭＣ乱数発生／基本統計
58
即ち、リープ・フロッグ法で新しく得られた変数
q
については、上に与えた確率 ( , | , ) q p q p
で採択の可否を決める。これは1に近い値のため、採択率はかなり高くなる。
最後に、オイラー法とリープ・フロッグ法の計算法を与えておく。
n
次元オイラー法( 1) ( )
( 1) ( ) ( )
p t p t h q
q t q t p t
n
次元リープ・フロッグ法(2 )
(2 2)
(2 1) (2 ) ( 2)
(2 2) (2 ) (2 1)
(2 2) (2 1) ( 2)
q t
q t
p t p t dh dq
q t p t p t
p t p t dh dq
実際の計算でのHMC法の使い勝手はどうであろうか。標準正規乱数を発生させるごとにリープ・
フロッグ法を用いるため、計算量（特に微分の部分）がかなり多くなる。採択率は上がるが、その分 計算量が増えるため、計算時間は MH 法に比べて長くなっている。しかし、乱数の精度から見ると 改良されているのではないかと思われる。
次に初期値が最尤値から離れている場合の収束性について、これまでの計算では、遠くの最尤値ま で速く収束するようには感じられない。むしろ初期値に対して最尤値が離れている場合は、密度関数 が計算誤差で0となってしまうところが問題のように思われる。これについてはMH法もHMC法 もあまり変わらないように思う。
59
図1 MCMC乱数発生メニュー
乱数発生には基本的にメトロポリス・ヘイスティングス法とハミルトニアン・モンテカルロ法があ り、メニューで変更できる。
プログラムを利用する際、まず「密度関数」テキストボックスに、出力させる目的分布の乱数の密 度関数を入力する。「例」のコンボボックスにサンプルが入っているので、それを参考にしてもらい たい。ここではまず、密度関数 = 1/6*exp(-abs(x)/3) の1次元の例を用いて説明を行う。
目的分布の密度関数を入力したら、描画範囲のx軸の上限と下限を入力する。この範囲はあくまで 描画する際の表示範囲で、乱数発生はこれにとらわれない。乱数の発生範囲は、「最大・最小」ボタ ンで、図2のように表示される。
図2 乱数発生の最小・最大 描画範囲が不明の場合はこの結果を参考にしてもよい。
描画範囲として下限-20と上限20を入力したら、まず、「ヒストグラム」ボタンで図3aのようなヒ ストグラムを描いてみる。
ＭＣＭＣ乱数発生／基本統計
60
図3a 乱数のヒストグラムと理論曲線（Seed=1）
ヒストグラムと同時に出力した乱数の統計量も表示される。採択率は、Metropolis-Hastingsアル ゴリズムの抽出率をいう。
図3aの中の曲線は目的分布の密度関数を利用した理論値である。この場合少しずれているが、乱 数の「Seed」を変えることによって分布が異なってくる。例として、図3bにSeed = 2 の場合を示 す。
図3b 乱数のヒストグラムと理論曲線（Seed=2）
ヒストグラムの階級幅は「x分割」の数によって決まる。この場合、範囲が40でx分割数が20で あるので階級幅は2になっている。
密度関数の形は、「描画」ボタンで見ることができる。但し、１変量関数グラフのプログラムを利 用するので、そのメニューが表示されるが、その中の「グラフ描画」ボタンをクリックすると図4の ようなグラフが表示される。
図4 密度関数グラフ