マルコフ連鎖モンテカルロ（MCMC）：複雑な確率分布からのサンプリング

確率分布の中には、直接扱うにはあまりに複雑なものがあります。

実世界のデータをモデリングすると、役に立つ結果に到達する前にしばしば数学が破綻します。積分は紙の上では扱えそうに見えても、潜在変数を少し追加した瞬間に手に負えなくなることがよくあります。とりわけベイズ推論では、事後分布が事前信念と観測データを組み合わせるため、単純な式で要約できない形になるのが一般的です。

マルコフ連鎖モンテカルロ（MCMC）の基本的な考え方は、数式を直接いじる代わりに、シミュレーションで分布を探索することです。分布全体を計算しなくても、その形状を反映したサンプルを生成します。

本記事では、MCMCの中核概念を解説し、代表的なアルゴリズムを一通り紹介し、Pythonでの実践方法を示します。

Python の数学をおさらいしますか？Demystifying Mathematical Concepts for Deep Learning のブログ記事を読み、Numpy での数学の適用例をご確認ください。

マルコフ連鎖モンテカルロとは？

マルコフ連鎖モンテカルロ（MCMC）は、確率分布が直接扱えないほど複雑な場合でも、その分布からサンプルを生成するためのアルゴリズム群です。

名前は2つの部分から成ります。マルコフ連鎖は、アルゴリズムが取り得る状態間をどのように遷移するかを制御します。各ステップは、それまでの履歴全体ではなく、今いる場所のみに依存します。モンテカルロは、関心のある量を推定するために乱数サンプリングを用いることを意味します。

両者を組み合わせることで、MCMCはランダムサンプルの連鎖を構築し、時間とともに目的の分布の形状を反映します。これは何よりもまずサンプリング手法です。数学を厳密に解くのではなく、シミュレーションで近似します。

MCMCが必要とされる理由

現実世界のデータ分布は、教科書的な分布ほどきれいではありません。

ベイズ推論では、観測データを踏まえたモデルパラメータの更新後確率、すなわち事後分布の計算をしばしば目指します。式は紙の上では簡単に見えます。事前分布に尤度を掛け、周辺尤度で割るだけです。最後の項は、あり得るすべてのパラメータ値にわたる積分を要します。高次元では、その積分は実質的に計算不可能です。

モデルが大きくなるほど状況は悪化します。パラメータや潜在変数を増やすと、厳密計算は袋小路になります。次のような一般的な場面でこれに直面します：

• ベイズ階層モデル： パラメータが他のパラメータに依存し、入れ子の分布となってきれいに積分できない
• 高次元の事後分布： パラメータ数が多すぎて格子探索の近似が計算的に現実的でない
• 非共役事前分布： 事前分布と尤度を組み合わせても既知の分布にならない

こうした場面でMCMCは有効な回避策です。分布を計算する代わりに、そこからサンプルを引きます。サンプルが、積分を解かずとも必要な情報をすべて含みます。

MCMCを支える2つのアイデア

MCMCは、それ自体は単純だが組み合わせると強力な2つの考え方を統合しています。順を追って説明します。

マルコフ連鎖

マルコフ連鎖は、各ステップが現在の状態のみに依存する状態の列です。

どこにいたかは重要ではありません。次にどこへ行くかは現在の状態だけで決まります。この「記憶なし」の性質（形式的にはマルコフ性）が、数学を扱いやすくし、アルゴリズムを実用的にします。

連鎖は1ステップずつ状態を移動し、適切に設定すれば最終的に定常分布に落ち着きます。これは、どの状態にいる確率も変化しなくなる安定したパターンです。MCMCはまさにこの安定分布を利用するよう設計されています。

モンテカルロ法

モンテカルロ法は、直接計算が難しい量を乱数サンプリングで推定します。

分布から十分な数のランダムサンプルを引き、そのサンプルから平均や分散などの特性を推定します。サンプル数が増えるほど、推定値は真値に近づきます。

モンテカルロ法単体では、分布から直接サンプリングできることが前提になります——それこそが今解決したい問題です。ここをマルコフ連鎖が担います。

マルコフ連鎖モンテカルロの動作原理

MCMCは、各ステップで単純な意思決定を行うループです。

1. 初期状態を決める： パラメータ空間の出発点を選びます。大まかな当て推量でも構いません
2. 新しい状態を提案する： 提案メカニズムを使って次の移動先を示します
3. 受け入れるかを判断する： 提案状態と現在状態を比較し、目的分布の下での各状態の尤もらしさに基づいて受容または棄却します
4. 反復： このループを何千回も回し、その過程でサンプルを収集します
5. サンプルを利用する： 収集したサンプルから平均や信用区間など関心のある量を推定します

受容・棄却のステップが「魔法」が起きる箇所です。

良い状態は悪い状態より高い頻度で受け入れることで、連鎖は高確率領域へ引き寄せられます——分布全体を計算する必要はありません。

初期のサンプルは開始点に依存するため破棄します。十分な反復後には、連鎖は開始点を忘れ、残ったサンプルが目的分布の真の形を反映します。

目的分布と定常分布

MCMCは、直接サンプリングできない目的分布からサンプルを生成することを目標に構成されています。

目的分布は、通常はベイズ推論における事後分布など、学びたい対象です。正規化定数を除いた形は分かっていても、その定数自体は直接計算できません。MCMCはそれを必要としません。

すべてのMCMCアルゴリズムは、自身のマルコフ連鎖が目的分布と一致する定常分布を持つよう設計されています。定常分布とは、十分なステップ後に連鎖が落ち着く分布のことです。

連鎖を走らせ続ければ、目的分布からの標本と見分けがつかないサンプルが得られます。積分は回避されます。

代表的なMCMCアルゴリズム

実務でよく目にするMCMCアルゴリズムはいくつかあります。どれも中核のループは同じですが、新しい状態の提案方法や目的分布の情報の使い方が異なります。

メトロポリス法

メトロポリス法は最も単純で、すべての始まりとなったMCMCアルゴリズムです。

各ステップで、現在の状態にランダムなノイズを加えて新しい状態を提案します。提案状態が目的分布の下でより高い確率なら常に受け入れ、低い場合は2つの確率の比に比例した確率で受け入れます。そうでなければ連鎖はそのまま留まります。

この受容・棄却メカニズムにより、分布全体を計算せずとも、連鎖は高確率領域に長く滞在します。

メトロポリス法は対称な提案分布を用い、どの方向にも同じ確率でステップを提案します。モデルが大きくなると破綻しがちです。

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

メトロポリス–ヘイスティングス（MH）法は、非対称な提案分布を許容することでメトロポリス法を一般化します。

MHは、ある提案が他よりも起こりやすい事実を受容確率に反映します。提案を目的分布の形に合わせて調整でき、探索が改善し、収束が速くなります。

現代のMCMC手法の多くは、MHの拡張か、同じ原理に基づいています。メトロポリス–ヘイスティングスを理解すれば、この分野の土台が分かります。

ギブスサンプリング

ギブスサンプリングは、すべてのパラメータを一度に更新する代わりに、1つずつ変数を更新します。

各ステップで、他の変数の現在値を条件とした条件付き分布から各変数をサンプリングします。すべての変数を一巡すると、1回の反復が完了です。

この方法では受容・棄却ステップが不要で、各条件付きサンプルは常に受け入れられます。結合分布からの直接サンプリングが難しくても、条件付き分布は扱える場合に有用で、ベイズ階層モデルで一般的です。

ハミルトニアン・モンテカルロ

ハミルトニアン・モンテカルロ（HMC）は、現代のベイズ推論をスケールさせた最初のアルゴリズムです。

HMCは、ランダムに新状態を提案する代わりに、目的分布の勾配情報を用いて、現在位置から遠いが受け入れられやすい状態を提案します。ランダムウォーク法よりもはるかに効率よくパラメータ空間を移動し、棄却が少なく、高次元分布の探索が改善します。

メトロポリスのようなランダムウォーク法は、パラメータ数が増えるとスケールしません。HMCは同程度の問題を抱えません。

HMCは、最も広く使われる確率的プログラミング基盤の1つであるStanのエンジンです。PyMCで使われるNo-U-Turn Sampler（NUTS）は、HMCの適応型拡張で、ステップサイズやステップ数の手動チューニングを不要にします。

ベイズ統計におけるMCMC

MCMCが最大の影響を与えた領域を1つ挙げるなら、それはベイズ推論です。

ベイズ統計の中心は事後分布、すなわちデータを見た後のモデルパラメータの更新後確率です。これを計算するには事前分布に尤度を掛けて正規化しますが、この正規化には、多くの場合扱えない積分が必要です。

MCMCはこのステップを完全に取り除きます。正規化されていない事後分布を任意の点で評価し、あとは連鎖に任せます。

簡単な例として、コインの偏りを推定するとしましょう。コインはおそらく公平だという事前信念から始め、コイントスの結果列を観測します。単純なコインモデルなら事後分布は閉じた形を持ちますが、100枚のコインの偏りを同時に推定するような階層構造を加えると、計算が不可能になります。

MCMCなら、連鎖を回して事後分布からサンプルを集め、それを使って必要な量を計算できます。

バーンイン、収束、ミキシングの理解

この3つの概念は、MCMC初心者のデータサイエンティストを混乱させます。ここを誤ると結果は得られても、なぜ信頼できないのか分からなくなります。

バーンイン（Burn-in）

マルコフ連鎖の開始時点では、目的分布の高確率領域がどこか分かりません。

初期サンプルは目的分布ではなく開始点に影響されます。バーンインとはそれらを破棄する慣行です。まず一定回数連鎖を走らせ、そのサンプルを捨て、連鎖が良い出発点を見つけた後のサンプルだけを残します。

バーンインの長さに普遍的な規則はありません。モデルや開始点、ミキシングの良し悪しに依存します。実務では、あらかじめ固定値を決めるよりも、トレースプロットで視覚的に診断します。

収束（Convergence）

収束とは、連鎖が開始点の影響を受けなくなり、目的分布を反映したサンプルを生成している状態です。

収束していない連鎖はバイアスのあるサンプルを生みます。そこから計算した平均は真の事後平均に一致せず、連鎖がたまたま滞留していた場所を反映します。

収束は診断により事後的に評価します。異なる初期値から複数の連鎖を走らせ、結果が一致するか確認するのが、収束失敗を検出する最も信頼できる方法の1つです。

ミキシング（Mixing）

収束していてもミキシングが悪ければ問題は残ります。

ミキシングは、連鎖が目的分布をどれだけよく探索できているかを表します。ミキシングが良い連鎖は自由に動き回り、高確率・低確率の両領域を訪れ、ほぼ独立したサンプルを生みます。ミキシングが悪い連鎖は、移動までに多くの反復を要して1つの領域に留まり、相関の高いサンプルを生み、分布全体を代表しません。

ミキシングの悪さは、ノイズの多い水平帯ではなく、ゆっくり蛇行する川のようなトレースプロットとして現れます。そう見えたらサンプラーのチューニングが必要です——より良い提案分布にするか、アルゴリズム自体を変えるべきです。

ミキシング比較プロット

MCMC結果の評価方法

ここでは、MCMCを評価する4つの方法と、それぞれを使う場面を説明します。

トレースプロット

トレースプロットは、各反復におけるパラメータのサンプル値を示します。MCMCを実行した後、最初に確認すべきものです。

健全なトレースプロットは、安定した平均の周りの白色雑音のように見えます。トレンド、長いフラットな区間、ゆっくりしたドリフトは見られないはずです。連鎖がさまよったり、1つの領域に多くの反復で留まっているなら、それはミキシングの問題であり、サンプルは信頼できません。

トレースプロットの可視化

自己相関

MCMCのサンプルは決して完全に独立ではありません。各サンプルは直前のサンプルの影響を受けます。自己相関は、反復にまたがるサンプル間の相関の強さを測ります。

自己相関が高いと、サンプル数が示唆するほど情報量がありません。2000個の相関したサンプルが、200個の独立サンプルと同程度の情報しか持たないこともあります。多くのMCMCライブラリには自己相関プロットがあり、サンプル間隔が広がるにつれて相関がどれだけ速く低下するかを確認できます。

自己相関プロットの可視化

有効サンプルサイズ

有効サンプルサイズ（ESS）は、自己相関を実用的な数値に翻訳します。すなわち、その連鎖がいくつの独立サンプルに相当するかです。

5000サンプルを引いてもESSが200なら、実質的には200個の独立標本に相当する統計的な力しかありません。ESSが低ければ、連鎖を長く回すか、サンプラーを調整するか、その両方が必要です。実務では、各パラメータで少なくとも数百のESSを確保してから推定を信頼するのが一般的です。

収束診断

複数の連鎖を走らせると、それらが同じ分布に収束したかを形式的に検定できます。Gelman–Rubin 診断（R-hatとして報告）は、各連鎖内の分散と連鎖間の分散を比較します。

R-hatが1.0に近ければ連鎖の一致を示し、良いサインです。1.01や1.05を超える値（ライブラリの閾値に依存）は、収束しておらず、さらに反復が必要であることを示唆します。PyMCのような最新ライブラリはR-hatを自動計算し、値が高すぎると警告を出します。

PythonでのMCMC

Pythonには、思想の異なるMCMCライブラリがいくつかあります。

• PyMC： 多くのデータサイエンティストにとって最も取り組みやすい選択肢です。Python構文でモデルを定義すると、PyMCが内部でNUTS（HMCの適応版）を使ってサンプリングを実行します。Pythonでのベイズモデリングの定番です。
• CmdStanPy / PyStan： 専用の確率的プログラミング言語であるStanへのPythonインターフェースです。Stanはモデルを実行前にC++へコンパイルするため高速で、複雑な階層モデルに適しています。その代わり学習曲線は急です。
• NumPy： 専用のMCMCライブラリはありませんが、メトロポリス–ヘイスティングスのようなアルゴリズムを一から実装できます。高水準ツールに進む前に、内部動作を理解するのに役立ちます。

多くの実務では、PyMCから始めるのが良いでしょう。ここでもそれを使用します。追従する場合は、まずライブラリをインストールしてください：

pip install pymc

話を簡単にするため、コイントスの列からコインの偏りを推定するという易しい例に絞ります。

ステップ1：モデルを定義する

import pymc as pm
import numpy as np

# 1 = heads, 0 = tails
observed_flips = np.array([1, 0, 1, 1, 0, 1, 1, 1, 0, 1])

with pm.Model() as coin_model:
    # Prior: we believe the coin is probably fair
    bias = pm.Beta("bias", alpha=2, beta=2)

    # Likelihood: observed flips given the bias
    flips = pm.Bernoulli("flips", p=bias, observed=observed_flips)

pm.Betaの事前分布は、コインが公平だという弱い信念を表現します。pm.Bernoulliの尤度は、モデルを観測データに結びつけます。

ステップ2：サンプラーを実行する

with coin_model:
    trace = pm.sample(2000, tune=1000, return_inferencedata=True)

サンプラー実行の出力

tuneはバーンインのステップ数を設定し——それらのサンプルは破棄されます。sampleは、チューニング後にチェーンごとに2000個の事後サンプルを引きます。

ステップ3：事後分布を確認する

import arviz as az

az.plot_trace(trace, var_names=["bias"])
az.summary(trace, var_names=["bias"])

モデルのトレースプロットと要約結果

az.summary()は、各パラメータの事後平均、標準偏差、R-hatを返します。R-hatが1.0に近ければ、連鎖は収束しています。az.plot_trace()は各パラメータについて、トレースと事後分布を並べて描画します。

このデータセット（10回中7回の表）では、事後平均は0.642、標準偏差は0.124です。これはデータの証拠を反映しつつ、公平なコインという事前に近い値になっています。R-hatは1.00、ESSは2000を大きく上回っており、連鎖は収束し、サンプルは信頼できます。

MCMCでよくあるミス

MCMCは実行自体は簡単ですが、誤用もしやすいです。よく見られるミスを挙げます。

• 早すぎる収束の思い込み：サンプラーがエラーなく走ったからといって、収束したとは限りません。結果を信頼する前に、必ずR-hatとトレースプロットを確認してください。見かけ上安定していても、パラメータ空間の一領域に引っかかっている場合があります。
• バーンインの無視：初期サンプルは目的分布ではなく開始点に偏っています。破棄してください。多くのライブラリはtuneパラメータで自動処理しますが、分析に誤って含めていないか確認しましょう。
• ミキシング不良： 動きが遅い連鎖は強く相関したサンプルを生み、数ほどの情報を運びません。トレースプロットがノイズの多い水平帯ではなく、ゆっくりしたドリフトに見えるなら、提案分布の調整か、より良いアルゴリズムが必要です。メトロポリスからNUTSに切り替えるだけで解決することも多いです。
• サンプル数が少なすぎる： 有効サンプルサイズ（ESS）が低いと、連鎖が形式上は収束していても推定は信頼できません。ESSが総ドロー数のごく一部しかないなら、連鎖を長く回すか、ミキシングを改善してから結論を出してください。

MCMCと他のサンプリング手法の比較

分布を近似する方法はMCMCだけではありません。代替手法との比較は次のとおりです。

• 棄却サンプリングは、より単純な分布から提案を引き、目的分布への適合度に基づいて受け入れます。うまく機能すれば厳密ですが、高次元では受容率が崩壊します。
• 重要サンプリングは、提案分布からのサンプルに重み付けを行い、目的分布下での期待値を近似します。低次元では高速かつ有効ですが、次元が増えると重みの分散が発散します。実際には、わずかなサンプルがほとんどの重みを担い、推定が不安定になります。
• 変分推論（VI）は、単純な分布を当てはめて目的分布を近似し、あるダイバージェンスを最小化します。VIはMCMCより高速で、大規模データにもスケールしますが、近似誤差を伴います。MCMCなら捉えられる事後構造を見逃す可能性があります。

要点は次のとおりです：

MCMCと代替手法の比較

精度が速度より重要なとき、MCMCが適切です。大規模データにスケールさせる必要がある、またはリアルタイム推論が必要な場合は、精度低下のトレードオフを受け入れて変分推論が有力になります。

まとめ

MCMCは外から見ると難しそうですが、実際にしていること——直接計算できない分布の形を徐々に反映するサンプルの連鎖を構築する——を理解すれば、納得のいくツールです。

また、マルコフ連鎖とモンテカルロ法に分解して考えると、いっそう理解しやすくなります。

ベイズ統計における役割は計り知れません。通常は手の届かない事後分布も、信頼できるサンプラーを得た瞬間に解けるようになります。だからこそ、MCMCはPyMCやStanのような確率的プログラミングライブラリの中核なのです。

ただし実装に飛び込む前に、直感を固めておきましょう。連鎖がバーンインを必要とする理由、ミキシングの意味、トレースプロットの読み方を理解してください。コード自体は簡単で、Pythonライブラリが抽象化をシンプルな関数呼び出しの背後に隠してくれます。

機械学習に熟達したい方は、Machine Learning Scientist in Python トラックにご登録ください。85時間の教材で、2026年の就業準備が整います。