確率的プログラミング
Probabilistic Programming ・ かくりつてきぷろぐらみんぐ
確率モデルをコードで記述し、事後分布の推論を処理系に任せるパラダイム
概要
確率的プログラミングは、確率モデルをプログラミング言語のコードとして記述し、事後分布の計算という難しい部分を処理系(推論エンジン)に任せるパラダイムです。そのための言語・ライブラリを確率的プログラミング言語(PPL: Probabilistic Programming Language)と呼び、Stan、PyMC、NumPyro、Turing.jl などが代表格です。
書き手がやることは「データがどう生成されたと考えるか」を宣言することだけです。パラメータに事前分布を置き、データがパラメータからどんな分布で生成されるかを書く — この数行のモデル記述は、統計の教科書に載る数式とほぼ一対一に対応します。実行すると処理系がMCMCなどのアルゴリズムを走らせ、事後分布からのサンプルを返してくれます。「モデルを考える仕事」と「推論アルゴリズムを実装する仕事」の分離こそが、このパラダイムの核心です。
なぜ生まれたか
ベイズ推論の実務には長らく高い壁がありました。共役事前分布で解析的に解ける教科書的なモデルを一歩でも外れると、事後分布は手計算では求まらず、モデルごとに専用のMCMCサンプラーを自分で導出・実装する必要があったのです。モデルを少し変えるたびに数週間のプログラミングをやり直す — 統計学者の時間の大半が、統計ではなく数値計算の実装に費やされていました。
これを最初に大きく変えたのが1990年代のBUGS(Bayesian inference Using Gibbs Sampling)です。モデルの構造を宣言的に書けば、ギブスサンプリングを自動で組み立ててくれる — 「モデルを書けば推論はついてくる」という体験を初めて広めました。しかしギブスサンプリングは階層ベイズモデルのようにパラメータ間の相関が強いモデルで効率が悪く、限界がありました。2012年に登場した Stan は、勾配情報を使うハミルトニアン・モンテカルロ(HMC)とその自動チューニング版 NUTS、そして勾配を機械的に計算する自動微分を組み合わせ、数百〜数千パラメータのモデルでも実用時間で推論できるようにしました。同時期以降、深層学習フレームワークの自動微分基盤に乗った PyMC や NumPyro が広がり、確率的プログラミングはベイズ統計の標準的な実践方法になっています。
詳細
モデル記述はどう見えるか
たとえば「コインを20回投げて13回表だった。表の出る確率 p はどれくらいか」というモデルは、PPL ではおおむね次の2行に対応します。「p はベータ分布 Beta(1,1) に従う」「表の回数は二項分布 Binomial(20, p) に従い、観測値は13」。前者が事前分布、後者が尤度で、これは生成モデル — データが作られる手順の記述 — をそのままコードにしたものです。実行すれば p の事後分布からのサンプルが数千個返り、平均を取れば点推定、分位点を取れば信用区間が得られます。「事後分布を数式で求める」代わりに「事後分布からのサンプルの集まりとして持つ」のが、この世界の基本的な作法です。
処理系の中で起きていること
ワークフロー全体は次のように流れます。
鍵になっている技術は2つです。1つは自動微分で、書かれたモデルから対数事後密度の勾配を機械的に導出します。もう1つはその勾配を使うHMC/NUTSサンプラーで、事後分布の「地形」を物理シミュレーションのように滑走して、酔歩型のマルコフ連鎖よりはるかに効率よく探索します。この2つの組み合わせにより、書き手は微分もサンプラー設計も意識せずに済みます。サンプルが得られたあとの期待値や確率の計算は、すべてモンテカルロ法の平均計算です。
「推論の自動化」がモデリングを変えた
推論が自動化されたことの意味は、単なる省力化にとどまりません。モデルの修正と再実行が数分で回るため、「単純なモデルで始め、データとのずれを見て少しずつ拡張する」という反復的なモデリングが現実的になりました。外れ値が多ければ正規分布を裾の重いt分布に差し替える、グループ構造があれば階層ベイズモデルにする、欠測があれば欠損データを未知パラメータとして一緒に推定する — いずれもモデル記述を数行変えるだけで、推論側の書き直しは不要です。頻度論の枠組みでは手法名(t検定、分散分析…)を選ぶ形になりがちな分析が、「生成過程を自分で書く」自由なモデリングに変わるのです。
使いどころと落とし穴
MCMCが標準の推論手段ですが、データが大きい・速度が要るときには変分推論(事後分布を扱いやすい分布で近似する最適化ベースの方法)を選べる処理系も多くあります。精度と速度のトレードオフで、変分推論は速い代わりに不確実性を過小評価しがちです。
落とし穴の筆頭は「動いた=正しい」ではないことです。サンプラーは必ず収束診断とセットで使います。複数のチェーンが同じ分布に収束したかを見る R-hat(1.01以下が目安)、実質的な情報量を測る有効サンプルサイズ、HMC特有の警告である divergence の有無を確認し、問題があればモデルの再パラメータ化や事前分布の見直しで対処します。また、推論が正しくてもモデルが現実に合っていなければ意味がないため、事後分布から擬似データを生成して実データと見比べる事後予測チェックを行うのが定石です。識別不能なパラメータ(データからは和しか決まらない2つの効果など)や無情報すぎる事前分布は、収束不良という形で表面化することが多く、診断はモデルの欠陥を教えてくれるセンサーでもあります。推論エンジンは計算を肩代わりしてくれますが、モデルの責任は書き手にある — これが確率的プログラミングと付き合う基本姿勢です。
