/w=3840,quality=90,fit=scale-down)
※本記事は、Stanford University提供の講義「Stanford AA228V: Validation of Safety Critical Systems」の内容を基に作成されています。講義の詳細情報は公式ウェブサイト https://aa228v.stanford.edu/ でご覧いただけます。また、講義で使用されている教科書は https://algorithmsbook.com/validation/ で入手可能です。本記事では、講義内容を要約しております。
この講義は、Sydney Katz氏（スタンフォード大学ポスドク研究員）によって提供されています。Sydney Katz氏についての詳細は https://sydneymkatz.com/ でご覧いただけます。
なお、本記事の内容は原著作者の見解を正確に反映するよう努めていますが、要約や解釈による誤りがある可能性もありますので、正確な情報や文脈については、オリジナルの講義動画をご視聴いただくことをお勧めいたします。この講義の受講登録や詳細については、Stanford Online（https://online.stanford.edu/ ）をご参照ください。Stanford Onlineはスタンフォード大学の各学部・部門が提供する学術的・専門的教育へのポータルであり、スタンフォード工学部のグローバル・オンライン教育センター（CGOE）によって運営・管理されています。
1. 導入
1.1. 講義の概要と期待
Sydney Katz: 今日の講義についてとても興奮しています。この講義の準備をしていたとき、今日のスライドの一部は私が最初に作ったものです。数学的にも視覚的にも非常に満足のいくトピックだと思うので、皆さんも楽しんでいただければと思います。
本講義を始める前に一言言っておきたいのは、これらのトピックは私が最初に学んだ時、本当に理解するまでに何度か繰り返して学ぶ必要がありました。もしかしたら、それは単に私がこれから皆さんに提示するほど美しく説明されていなかったからかもしれません。そのため、皆さんにとってはそれほど難しくないかもしれませんが、これらの概念は少し時間をかけて理解する必要があるかもしれません。この講義を見て、教科書に戻り、また講義を見直すことで、徐々に理解が深まるでしょう。
1.2. 失敗分布の重要性
Sydney Katz: この1週間ほど、私たちは「falsification」（偽証）というタスクについて話してきました。このタスクでは、システムの失敗を見つけることに焦点を当てています。場合によっては最も起こりやすい失敗を探したり、見つけられる失敗を何でも探したりしています。しかし、時には一つの最も可能性の高い失敗だけでなく、観測される可能性のある失敗の完全な分布を知りたい場合もあります。
今日は、この考え方をさらに一歩進めて、失敗の完全な分布からサンプリングする方法について話します。おそらく気づかれると思いますが、そのプロセスはやや複雑で難しいかもしれませんが、失敗について起こりうる分布全体を把握できるという追加の報酬が得られます。
2. 失敗分布の数学的定義
2.1. 軌道分布の定義
Sydney Katz: 数学的に「失敗分布」が実際に何であるかについて話しましょう。前回の講義で「軌道分布（trajectory distribution）」について触れました。これは特定の軌道が与えられたときに、その軌道がどれくらい可能性があるかを示すものです。
基本的に、この分布は観測される可能性のあるすべての軌道に対する分布を表しています。これを「p(τ)」と表記し、特定の軌道τの確率または可能性を示します。これは私たちが「名目軌道分布（nominal trajectory distribution）」と呼ぶもので、システムが通常どのように振る舞うかを表現しています。
2.2. 失敗分布の条件付き確率としての表現
Sydney Katz: 失敗分布は、この軌道分布を条件付き確率として拡張したものです。具体的には、失敗分布はτが失敗であるという条件の下での確率と定義します。数式で表すと「p(τ|τ∉ψ)」となります。ここで「τ∉ψ」という表記は、軌道τが仕様を満たさない、つまり失敗であることを意味します。
この分布は実際には次のように展開できます。まず分子に注目しましょう。中括弧で囲まれた部分は、このコースで使用している指示関数（indicator function）の表記です。この関数は、τが仕様を満たさない（つまり失敗である）場合に1を返し、成功である場合には0を返します。そしてこの「p(τ)」は名目軌道分布です。
つまり、この分子は、軌道が失敗の場合はp(τ)を返し、そうでない場合は0を返します。これが失敗分布の本質です。しかし、もう一つ追加する必要があるのが分母の部分です。この分母がないと、いくつかの軌道に0の確率を割り当てることになるため、確率が1に積分されなくなります。そのため、有効な確率分布を作るためには、この積分で割る必要があります。
つまり、すべての可能な軌道にわたってこの積分を行い、失敗分布を得ます。これを視覚的に見てみましょう。
2.3. 非正規化確率密度関数の概念
Sydney Katz: この失敗分布の式を見ると、一般的には分母部分、つまり積分部分の計算が非常に難しいことがわかります。場合によっては計算が不可能なこともあります。これは残念なことですが、通常、この失敗分布の解析的な表現を計算することはできません。
しかし、朗報があります。分子部分は計算できるのです。これは以前のベイズ推定について話した時にも似たような考え方がありました。特定の軌道が与えられたとき、この式の分子部分の計算方法を知っています。最初の項を計算するには、その軌道が失敗かどうかを確認するだけでよく、失敗であれば1を返し、そうでなければ0を返します。そして2番目の項を計算するには、名目軌道分布の下でその軌道の可能性を評価するだけです。前章で尤度目標などについて話したとき、特定の軌道が発生する確率の計算方法についてすでに説明しました。
つまり、この分子は計算可能なのです。前のスライドで示したように、分子は実際に分布の「形状」を決めるものです。したがって、この分子だけを知ることでも、興味深いことが行えるようになります。
これは「非正規化確率密度（unnormalized probability density）」と呼ばれます。これは確率密度の正しい形状を与えますが、1に積分されないため有効な確率分布ではありません。しかし、それでも情報を提供することができます。このクラスでは、非正規化確率密度について話すとき、変数の上にバーを付けて表記します。
この非正規化確率密度から実際にサンプルを取得できるアルゴリズムが存在するという点が非常に魅力的です。以前ベイズ推定について話した時に触れましたが、ある種の「魔法」が起こり、非正規化確率密度からサンプルを取得できるのです。今日はその「魔法」の正体について少し説明します。
このようにして、失敗分布を明示的に記述できなくても、その分布からサンプルを多数取得することで暗黙的に表現することができます。例えば、グリッドワールドの場合、その失敗分布からサンプリングして、様々な異なる軌道を取得し、それらをすべて一緒にプロットできます。これにより、グリッドワールドシステムの失敗分布がどのようなものかがわかります。
3. 1次元ガウス分布の例
3.1. 単純なシステムでの失敗分布の視覚化
Sydney Katz: これを視覚的に理解するために、非常に単純な例を使ってみましょう。プロジェクト1を始めた方は、これがプロジェクト1の小さなシステムで使用している例だと気づくでしょう。
これは軌道が基本的に1つの状態だけで構成される超シンプルなシステムです。初期状態をサンプリングするだけで、それが軌道のすべてです。この初期状態、つまり私たちのτは、このようなガウス分布からサンプリングされます。これが私たちの名目軌道分布「p(τ)」であり、このホワイトの分布で表されています。
この式の構成要素が、この分布に関してどのように見えるかを見てみましょう。ここでは、サンプリングした状態が特定の値（この場合はマイナス1）より小さい場合、その軌道は失敗であると定義します。つまり、左側にあるものはすべて失敗で、右側にあるものは失敗ではありません。
この式の分子部分がこの特定のシステムに関してどのように見えるかを見てみましょう。先ほど説明したように、この分子は軌道が失敗である場合、名目軌道分布と同じ値を持ちます。この赤い線は分子を示しており、失敗である軌道（マイナス1より小さいもの）に対しては名目軌道分布と全く同じ値を持ち、他のすべての軌道に対しては0の値を持ちます。
分母については、これはこの値のすべての可能な軌道にわたる積分です。微積分学から、この積分はこの曲線の下の面積であることがわかります。これは定数値なので、この分布を取り、この定数値で割ると考えることができます。この面積はτに基づいて変化するわけではなく、常に同じです。なぜなら、すべての可能なτにわたって計算しているからです。結果として、この値をスケーリングするだけで、曲線下の面積が1に等しくなります。これで有効な確率分布ができました。
つまり、ここで言いたいのは、基本的にこの分子が失敗分布の実際の形状を与えるということです。そしてこの分母は、その形状をスケーリングして1に積分される分布を得るために計算するものです。
3.2. 分布の形状と正規化定数の関係
Sydney Katz: 先ほどお見せしたように、上部の分子が分布の形状を与え、下部の分母が正規化定数と呼ばれるスケーリング係数を与えます。これは単なる定数で、曲線下の面積を与えるものです。
もしこの計算ができれば、すべてのシステムに対する失敗分布を知ることができます。しかし、ここでの主な問題は、通常、この分母、つまり積分部分の計算が非常に難しいということです。場合によっては計算が不可能なこともあります。つまり、失敗分布の解析的な表現を計算できないことが多いのです。
しかし、朗報があります。分子部分は計算できます。クラスの前半でベイズ推定について話したときに似たような考え方がありました。特定の軌道が与えられたとき、この式の分子の計算方法は知っています。最初の項を計算するには、その軌道が失敗かどうかをチェックするだけでよく、失敗であれば1を返し、そうでなければ0を返します。そして2番目の項を計算するには、名目軌道分布の下でその軌道の可能性を評価するだけです。
つまり、この分子は計算可能です。そして前のスライドで示したように、分子は分布の形状を決めるものです。したがって、分子だけを知ることでも、興味深いことができます。これを「非正規化確率密度」と呼びます。これは確率密度に正しい形状を与えますが、1に積分されないため有効な確率分布ではありません。
実は、非正規化確率密度からサンプルを取り出すことができるアルゴリズムが存在します。これにより、失敗分布を明示的に書き出すことはできなくても、その分布からサンプルを多数取得することで暗黙的に表現できるのです。例えば、グリッドワールドでは、その失敗分布からサンプリングし、様々な異なる軌道を取得し、それらをすべて一緒にプロットすることで、グリッドワールドシステムの失敗分布がどのようなものかを示すことができます。
このような1次元ガウスの例は単純ですが、これらの手法は実際にはこの1次元ガウスよりも大きなシステム、例えば50次元、100次元の振り子システムにも適用できます。これらの手法がどのようにスケールアップするかを見るのは非常に興味深いです。まずは1次元ガウスで直感を得て、後でより大きなシステムに移行します。
4. リジェクションサンプリング
4.1. ダーツボードアナロジーによる直感的説明
Sydney Katz: それでは、リジェクションサンプリングについて説明しましょう。リジェクションサンプリングを理解するために、まずダーツボードのアナロジーを使います。
長方形のダーツボードがあると想像してください。通常のダーツボードとは少し違い、単なる大きな長方形で、そこにダーツを投げるとボードに刺さります。
最初にやることは、このダーツボード上にサンプリングしたい目標の非正規化密度の形状を描くことです。これは先ほど説明した分子の部分で計算できる密度の形状です。そしてボード上に描きます。
次に、このボードに対して一様にダーツを投げます。つまり、何らかの方法で実際にランダムにダーツを投げ、ボード上のランダムな場所に着地させることができると仮定します。
そして、投げたすべてのダーツを見て、ボード上に描いた目標密度の上に着地したダーツと、その密度の下に着地したダーツを分離します。ここで「拒否」の部分が登場します。目標密度の上に着地したダーツをすべて拒否します。
そして、残りのダーツ（密度の下に着地したもの）をボードの下部に真っ直ぐに落とします。すると、これらのサンプルは実際にこの密度に従って分布しているのです。これは驚くべきことです。
これを視覚的に理解するのは難しいかもしれませんので、ビン（区間）を作成して、各ビンに落ちたサンプルをすべて積み上げてみましょう。そうすると、これらが実際に目標密度に従って分布していることがわかります。さらに、投げるダーツの数を増やせば、より良い近似が得られます。
この手法は通常の分布だけでなく、任意の形状の密度からサンプリングするために使用できます。例えば、別の形状の密度があったとしても、同じ手順でダーツをボードに投げ、密度の上下にあるものを分離し、上にあるものを拒否し、すべてを下に落とします。これらは実際にその特定の分布に従って分布しています。これをビンにまとめると、その関係がより明確に見えます。
4.2. リジェクションサンプリングの数学的定式化
Sydney Katz: これまでの直感的な説明をより数学的に定式化してみましょう。このアルゴリズムをより一般的に使えるよう、数学的表記を徐々に導入していきます。
通常の分布に戻りましょう。理解しやすいからです。まずダーツボードに軸を定義します。ダーツボードの横軸を変数τとします。これは軌道、つまり実際に抽出したいサンプルを表します。例えば、この軸が-3から3までの範囲だとします。そして縦軸をHとしましょう。これはボード上の高さを表します。範囲は0から0.6とします。
では、この分布からサンプルを抽出するために行った手順を数学的に書き出してみましょう。最初に行ったのはボードに一様にダーツを投げることでした。実際にダーツボードを使う場合、ダーツを投げると、非常に優れた技術を持っている人なら一様に投げることができるかもしれません。しかし、ここでは2段階のプロセスに分解します。最初は少し奇妙に思えるかもしれませんが、数学的に記述できる一般的なアルゴリズムに向けて進んでいるのです。
まず行うのは、ダーツボードのこの軸に沿って値をサンプリングすることです。-3から3の間で一様にサンプリングします。次に、値を選んだら、ダーツボード上で一様に高さを選びます。つまり、0から0.6の間で一様に値をサンプリングして、ダーツの高さを得ます。
この過程を多数のダーツに対して繰り返せば—例えば2番目のダーツでも同じことをして、-3から3の間で値をサンプリングし、再び0から0.6の間で高さをサンプリングします—これはボードに一様にダーツを投げるのと同じプロセスになります。2つのステップに分解しただけですが、各ダーツの位置を得るには同等のプロセスです。これを何度も繰り返すと、ボード全体に一様に分布したダーツが得られます。
そして、青い密度の上にあるダーツをすべて拒否し、残りのダーツをボードの下部に移動させます。このプロセスをもう一度行い、もう少し数学を加えます。
今度は、サンプリングしようとしている青い密度に名前を付けます。p̄(τ)と呼びましょう。非正規化密度を表すときはこのバーを使います。最初にやるのは-3から3の間で一様にサンプリングすることです。これを数学的に表現するために、実際にサンプリングしている分布を指定します。τはq(τ)という分布からサンプルされると言います。この場合、ボードの下部に沿って一様にサンプリングしているので、q(τ)は-3から3の間の一様分布です。これは先ほどと同じプロセスですが、サンプルが実際に由来する分布を定義しているだけです。
次に行ったのは、ボードの下部から高さまでの間で一様にサンプリングすることでした。このプロセスを記述する方法は、最初は少し奇妙に見えるかもしれません。ここで言っているのは、hが定数cとq(τ)の積に、0から1の間の一様分布からサンプリングした値rを掛けたものに等しくなるということです。定数cは、q(τ)をボードの高さまで引き上げるように選びます。つまり、c×q(τ)がダーツボードの高さを与えます。
これに0から1の間の一様分布からサンプリングした値rを掛けると、ダーツのy位置またはh位置が得られます。これは同じプロセスです、いくつかの変数を追加して何が起きているかを説明しただけです。
次に、このプロセスを使ってすべてのダーツをサンプリングしたら、青い密度よりも上にあるダーツを拒否します。数学的には、サンプリングした高さの値が、サンプリングしようとしている非正規化目標密度よりも高いかどうかをチェックするだけです。そうであれば拒否します。
最後に、残りのダーツをボードの下部に移動させます。そのためには、軸上の値だけが必要で、高さはもう気にしません。つまり、サンプルのτ値だけを保持します。これでターゲット分布からのサンプルが得られました。
ここで注意したいのは、これらの灰色のサンプルをすべて拒否し、ターゲット分布からのサンプルとして残ったのは青いサンプルだけだということです。もしターゲット分布から100個のサンプルが欲しいとしたら、100個以上のダーツを投げる必要があります。なぜなら、いくつかを拒否するからです。すべての灰色のダーツは無駄なサンプルと考えることができます。これらは結局、ターゲット分布のサンプルには寄与しませんでした。
これは、多くのサンプルを無駄にする場合、非効率になる可能性があります。もし代わりにこのようなダーツボードからサンプリングできたら、無駄になるサンプルがはるかに少なくなります。拒否されるサンプルが少なくなります。しかし、このダーツボードからダーツを確実にサンプリングできるようにする必要があります。重要なのは、これが可能だということです。実際、これを試みるべきです。常に完全に長方形のダーツボードのような一様分布を使いたいとは限りません。ただし、q(τ)をサンプリングできる分布にする必要があります。つまり、この紫色の分布からサンプリングできる必要があります。そして、その分布が常に目標より高くなるような値cを選ぶ必要があります。ターゲットより低い場合、必要な領域の一部を見逃してしまいます。
4.3. 提案分布の選択と効率性の関係
Sydney Katz: ここでダーツボードのアナロジーから離れて、より一般的にリジェクションサンプリングがどのようなものかについて話しましょう。アルゴリズムをもう一度書き直します。
アルゴリズムの最初のステップは、q(τ)から何らかのτをサンプリングすることです。ここでq(τ)は「提案分布（proposal distribution）」と呼ばれることが多いです。この「提案分布」という考え方は、今後数週間、失敗確率の推定について話す際に何度も出てきます。基本的に、提案分布とは、サンプリング方法が分からない分布からサンプルを得るために、サンプリング方法を知っている分布からサンプリングすることです。
そして基本的には、このダーツが特定の条件を満たす場合に拒否すると言いました。c×q(τ)×rがp̄より大きければ拒否します。これを少し書き換えると、両辺をc×q(τ)で割り、正にするために反転すると、rが特定の値より小さい場合にダーツを受け入れるということになります。
これが使用する基準です。rは単に0と1の間で一様にサンプリングした値であり、ほとんどのプログラミング言語ではrand関数を呼び出すだけでこれを行うことができます。
また、先ほど言ったように、ダーツボードが基本的に私たちの密度全体をカバーするのに十分な高さであることを確認する必要があります。つまり、提案分布q(τ)と値cを選択する必要があります。これは、サンプリングしようとしているターゲット密度が常にc×q(τ)以下であることを意味します。
アルゴリズムを簡単にまとめると、提案分布からサンプリングし、この確率に比例して受け入れます。そして、この条件を満たす分布を選ぶ必要があります。
さて、私たちの目標である失敗分布に戻りましょう。リジェクションサンプリングをこのシナリオに適用するには、ターゲット密度p̄(τ)が、τが仕様を満たさない（失敗である）という条件付きの密度である必要があります。そしてこれは、指示関数と名目軌道分布の積に等しいと言いました。
ここで、q(τ)、つまり提案分布を選択する必要があります。名目軌道分布を試してみましょう。つまり、q(τ)が名目軌道分布p(τ)に等しいとします。名目軌道分布からのサンプルの抽出方法は知っています。単に多数のロールアウトを実行するだけです。
そして値cを選択する必要があります。これが成り立つ必要があります。そこで代入してみましょう。p̄(τ)はこれに等しいとわかっています。したがって、この指示関数とp(τ)の積がc×q(τ)以下でなければなりません。q(τ)がp(τ)に等しいと言っているので、c×p(τ)です。
ここで、この不等式が真となるようなcの値を選ぶ必要があります。この場合、cとして1が適切です。1を代入すると、この指示関数とp(τ)の積がp(τ)以下であることを確認する必要があります。この指示関数は1か0なので、p(τ)か0のいずれかになります。p(τ)と0は常にp(τ)以下です。
したがって、c=1とすると、アルゴリズムはこのように簡略化されます：まずq(τ)からτをサンプリングします。q(τ)は名目軌道分布なので、名目からτをサンプリングします。そして、rがこのp̄/(c×q)より小さい場合にサンプルを受け入れます。これは指示関数と等しくなります。
これは少し複雑で奇妙に思えるかもしれませんが、実は非常に明白な結果を導きました。この指示関数は、τが失敗であれば1、そうでなければ0です。したがって、軌道を受け入れる確率は、それが失敗であれば1、そうでなければ0です。
つまり、名目分布からサンプリングし、失敗を得た場合はそれを保持し、成功を得た場合は捨てるということです。これは前の講義で、失敗分布からサンプリングするためのベースライン手法として紹介しました。
このすべてが複雑に見えましたが、実は非常に直感的な結果を導き出しました。名目分布からサンプリングし、失敗を保持するということです。そしてそれらの失敗は実際に失敗分布から来ることになります。
しかし、この方法には問題があります。名目分布を使った場合、失敗が稀なシステムでは非効率になる可能性があります。サンプルの多くを拒否し、あまり多くのサンプルが残らないからです。解決策の一つは、実際にサンプリングしようとしている分布に近い、より良い提案分布を考案することです。
4.4. 希少な失敗イベントに対する課題
Sydney Katz: 提案分布を名目分布にした場合の考え方を見てみましょう。これは直接サンプリングのアイデアと同じで、失敗でないものを拒否するだけです。サンプル数を増やすと、以前に示したダーツボードのような形になります。ダーツボードは基本的にこのような形をしています。このプロットは、真の失敗分布（計算可能なこの小さなおもちゃの問題のため）と、取得しているサンプルのヒストグラムを示しています。
サンプル数を増やしてより多くのサンプルを生成すると、その分布の形に大体近いものが得られます。しかし先ほど言ったように、これは稀な失敗イベントに対しては非常に非効率的かもしれません。閾値を左に移動させて失敗をより稀にしてみましょう。閾値を左に動かすと、この方法が崩れ始める点に達します。これらのサンプルはもはやこの密度に従って分布しているようには見えません。
より良い結果を得るために、手動で設計した新しい提案分布を選択した方が良いかもしれません。ここでは、この提案分布の平均を制御し、分布の高さを決める値cも制御します。平均をどの方向に動かすべきでしょうか？左に動かすべきです。（すでに0にいるので右には動かせませんでしたが）今は失敗閾値の周りにこの平均を中心に置きます。これが最も多くのサンプルを得られるように見えます。これですでに結果が改善し始めています。
このままでは非常に非効率的です。これは失敗密度をはるかに上回っています。もっと下に移動させた方が良いでしょう。cを減らすと、さらに良い結果が得られます。より少ないサンプルを拒否するため、分布を表現するためのサンプルがより多く得られます。
このcをこれ以上減らすとどうなると思いますか？白い部分より上に行くサンプルはカウントされなくなります。つまり、最大値の部分が過小評価されることになります。全プロセスが崩壊します。それは、十分に高くなければならないという要件が破られ、正しい分布が得られなくなります。
私たちがやりたいのは、cを可能な限り小さくすることですが、必要な特性を侵害しないようにします。しかし、一般的には、cが適切かどうかを判断するのは容易ではありません。超保守的な値を選ぶことがありますが、それでは役に立たないほど保守的になってしまう可能性があります。
名目軌道分布を使用すると、稀な失敗イベントを持つシステムに対して非効率になる可能性があることを見てきました。より良い提案分布を考案しようとすると、このプロセスを改善できますが、それを手動で設計するのは難しいかもしれません。これがリジェクションサンプリングの課題や欠点です。適切な提案分布を選択することは常に明白ではありません。また、適切な値cを選択することも簡単ではありません。
より大きな例（振り子など）でこれがどのように機能し、崩壊するかを見たい場合は、教科書の例6.1を参照してください。
5. マルコフ連鎖モンテカルロ法（MCMC）
5.1. MCMCの基本概念と直感的理解
Sydney Katz: 次のトピックであるマルコフ連鎖モンテカルロ法（MCMC）も同様に興味深いものです。これらの方法には並行性が見られます。物事を受け入れたり拒否したりします。しかし、これらの並行性を強く引き出そうとすると混乱するだけなので注意してください。目標は同じで、非正規化ターゲット密度からサンプリングすることですが、受け入れと拒否の方法は実際には同じではありません。
マルコフ連鎖モンテカルロについて話す前に、このトピックについて一言言っておきたいと思います。このトピックを初めて聞いた時、正直なところ少し恐れていました。なぜか怖いトピックのように思えました。しかし、これは非常に理解しやすいものだと言いたいです。特に実践的な観点からは。今日の目標は、このアルゴリズムの実践的な理解を提供することです。どのように使用し、何をするのか。これは非常に強力なものです。
そうは言っても、背後には深い統計理論的な意味や証明があります。それはこのクラスの範囲を超えています。しかし、この実践的な観点が、そちらに興味があれば深掘りするのに十分であることを願っています。そうでなければ、なぜ、どのように機能するかを知った上で、仕事で使用できます。
また、もし何らかのインポスターシンドロームを感じている人がいれば、「Markov」という言葉に最初に触れた私の経験を共有したいと思います。MITリンカーン研究所でのインターンシップの面接中、グループリーダーから私が関わる予定の仕事について説明がありました。そのグループはマイケルとロバートも所属していた、航空機の衝突回避に取り組むグループでした。彼らは特に「マルコフ決定過程」を使用していました。
面接中、後で覚えておけるよう、彼が言ったことすべてをメモしていました。最近、アパートを引っ越した際にそのメモを見つけましたが、昨夜はもう一度見せるために見つけることができませんでした。しかし、どんな感じだったか再現しました。「サーベイランスシステム、衝突回避…」などといった内容でしたが、メモには面白いものがありました。実は、そのメールの前にマルコフ決定過程について一度も聞いたことがありませんでした。そこで「mark off」と書いたのが私の試みでした。明らかに何が起きているのか分かっていませんでした。
後になって、マルコフは実際には多くのものに名前が付けられている有名な数学者であることを知りました。「mark off」という言葉ではありませんでした。もし少し怖いと感じているなら、私はここまで来ました。あなたもできます。あなたはすでにマルコフという言葉が何かを知っているので、より良いスタートを切っています。幸い彼らはこのメモを見なかったと思います。もし見ていたら、私を雇ったかどうか分かりません。しかし結果は良かったです。
マルコフ連鎖モンテカルロに入りましょう。マルコフ連鎖モンテカルロでは、サンプルを互いに独立して1つずつ抽出するのではなく、実際にサンプルのチェーンを維持します。マルコフ連鎖の考え方は基本的に、チェーン内の次のサンプルが前のサンプルに依存するチェーンです。
私たちが行うのは、何らかの初期サンプルから始めることです。何かを選ぶだけです。そして、次のサンプルが前のサンプルに依存するチェーンを構築します（どのような依存関係かはまだ説明していません）。それを繰り返して、この非常に長いチェーンを得ます。チェーンができたら、前と同様のことができます。すべてのサンプルを取得すると、このチェーン内のサンプルをすべて下に落とすと、サンプリングしようとしていた分布に従って分布することになります。
どのようにそれが機能するかについてはまだ説明していませんが、それが私たちが向かっているところです。基本的には、このサンプルのチェーンをどのように作成すれば、完成したときにこの分布に従って分布するかを理解したいのです。
5.2. サンプルチェーンの構築方法
Sydney Katz: では、このチェーンを実際にどのように作成するのでしょうか？まず、何らかの初期サンプルτから始める必要があります。何から始めるかは少し重要ですが、今はそれについて心配しないでおきましょう。つまり、何らかの初期軌道、初期サンプルから始めるということです。
アルゴリズムの各反復で、τに依存するこの条件付き分布から新しい「提案された」サンプルをサンプリングします。これを「τ'」（タウプライム）と呼びましょう。この分布はしばしば「カーネル」と呼ばれます。
具体的にするために、非常に一般的なカーネルは、チェーン内の現在のサンプルを中心としたガウス分布のようなものです。つまり、こんな感じです：現在のサンプルがあり、現在のサンプルを中心とした分布から新しいサンプルを引きます。例えば、このサンプルを引き、これが私たちのτ'です。
次に行うのは、リジェクションサンプリングと似ていますが、まったく同じというわけではありません。ある確率でこのサンプルを受け入れます。その確率は特にここに示す確率で、サンプリングしようとしている名目ターゲット密度とこのカーネル関数に依存します。この確率に関する直感をすぐに説明しますが、まずアルゴリズムの概要を把握しましょう。この確率に等しい確率で受け入れます。
つまり、randを呼び出すなどして、受け入れるかどうかを決定し、その確率に従って受け入れます。このサンプルを受け入れるとしましょう。それが次のチェーンのサンプルになります。
そして、カーネルに戻ります。これは再び、現在のサンプルを中心としたガウス分布のようなものです。新しい提案サンプルをサンプリングし、ここに示す確率に従ってそれを受け入れるかどうかを決定します。今回は受け入れないとしましょう。そうすると、現在持っているサンプルをそのまま保持します。
このプロセスを繰り返します。例えば、このサンプルを受け入れるなど、大きなチェーンを構築します。驚くべきことに、それが完了すると、そのチェーン内のサンプルは実際にサンプリングしようとしていたターゲット密度に従って分布しています。
素晴らしいことに、このプロセスは無限のサンプルの極限でターゲット分布に確実に収束します。永遠にサンプリングすれば、サンプリングしようとしていたターゲット分布と正確に一致します。
しかし、この魔法はどのように機能しているのでしょうか？なぜこのようになっているのでしょうか？直感を与えようと思います。それほど狂ったことではないのです。
まず最初に言いたいのは、カーネルは多くの場合、対称的に選ばれるということです。つまり、g(τ|τ')=g(τ'|τ)です。これは、これら2つが等しいことを意味します。だから、それらをキャンセルすることができます。これは奇妙で狂った仮定のように思えるかもしれませんが、ガウスカーネルが対称的であることを示す練習問題があります。だから奇妙なことではありません。ほとんどのカーネルは対称的なので、これを行うことができます。
しかし、これにより受け入れ確率を簡略化できます。今、私たちの受け入れ確率は、ターゲット密度の下での新しいサンプルの確率密度を、現在のサンプルの確率密度で割ったものになります。
すべてをわきに置いて考えてみてください。密度関数を渡され、「こんな感じです」と言われ、そこからサンプルを抽出するように言われたとしたら、その密度の高い領域からより頻繁にサンプルを引き、低い領域からは時々サンプルを引くでしょう。この受け入れ基準はそのような行動を強制します。
5.3. 受容確率の設定と根拠
Sydney Katz: ここには2つのケースがあります。最初のケースは、提案された新しいサンプルでのターゲット密度が、現在のサンプルでのターゲット密度よりも大きい場合です。この場合、τ'はτよりもターゲット密度の下でより可能性が高いと言っています。
これによりこの分数は1以上になります。つまり、1以上の確率で受け入れるということです。言い換えれば、常にそのサンプルを受け入れるということです。現在のものよりも高い可能性を持つ新しい提案サンプルを得るたびに、常にそれを受け入れます。
もう一つのケースは、新しいサンプルが現在のものよりもターゲット密度が低い場合です。これは「高密度からより頻繁にサンプリングする」という考え方を実現しています。もし低い場合、高密度からより頻繁にサンプリングしたいので、必ずしもそれを保持したくはありません。しかし、時々は低密度のサンプルも入れたいです。
基本的に言っているのは、τ'はターゲット密度の下でより可能性が低いので、1より低い確率でそれをサンプリングまたは受け入れるということです。その確率は、どれだけ可能性が低いかに比例します。もしそれがはるかに可能性が低ければ、それを受け入れる可能性はまだありますが、そのサンプルを受け入れる可能性は低くなります。
この方法が望ましいのは、高密度エリアからより頻繁にサンプリングするが、低密度エリアからも時々サンプルを得ることです。これは決してこれが機能するという証明ではありません。それはこのクラスの範囲を超えています。しかし、ここで起こっている魔法の一部を解明することを願っています。本当に起こっていることはそれだけです。
[質問：どうしてカーネルを対称にするのですか？]
カーネルについて、私たちはオールドサンプルの値に基づいて、その周りに確率分布を作成します。そして、それに基づいて新しいサンプルを得ます。しかし、ターゲットはオールドサンプルの周りの確率分布ではありません。
オールドサンプル周りの確率分布は、単にチェーンに対する新しい提案サンプルを得る方法です。以前のサンプルに条件付けられた分布からの抽出方法が必要です。それがこのカーネルの目的です。カーネルはpではなく、gです。
5.4. 対称カーネルの利点
Sydney Katz: 対称カーネルの主な利点は、受け入れ確率の式を大幅に簡略化できることです。カーネルが対称的であるとき、つまりg(τ|τ')=g(τ'|τ)であるとき、これらの項は相殺され、受け入れ確率は単にターゲット密度の比に簡略化されます。
これが対称カーネルを使用する理由の一つですが、一般的なガウスカーネルのような多くのカーネルは自然に対称的であるため、特別な設計は必要ありません。対称性があると、アルゴリズムがより理解しやすく、実装も簡単になります。
非対称カーネルを使用することもできますが、その場合は受け入れ確率の式がもっと複雑になります。gの項が相殺されないため、それらを計算に含める必要があります。
実際にMCMCを動作させるために、Plutoノートブックで例を見てみましょう。このプロットでは、チェーン内に一つのサンプルだけがある状態から始めています。これが選択した初期サンプルです。
このプロットはチェーンが時間とともにどのように構築されるかを示します。また、失敗領域も示しています。ここからサンプリングしようとしています。
2番目のプロットは現在のサンプルと提案されたサンプルを示しています。ガウスカーネルを使用していて、提案された次のサンプルは実際に確率が低く、分布からこのように遠くにサンプルを得ました。そして、これはサンプリングしようとしているターゲット密度を同じプロット上に示しています。
下のプロットでは、サンプルを追加するにつれて、見たすべてのサンプルのヒストグラムをターゲット密度と比較して示します。
次のサンプルに進む前に、この新しい提案サンプルを受け入れるか拒否するかを誰か教えてくれますか？
[回答：拒否します]
はい、拒否します。その理由は、新しいサンプルのターゲット密度が実際に0だからです。ここからサンプリングしようとしていて、もしより可能性が高ければ確率1で受け入れますが、明らかにより可能性が高くはありません。それは0です。もしより可能性が低ければ、新しい可能性と古い可能性の比率でサンプリングします。新しい可能性は0なので、このサンプルを受け入れる確率は0です。
次の反復に進むと、同じサンプルに留まっていることがわかります。今回は、このサンプルを受け入れる可能性があります。このバーは実際にその可能性がどれくらいであるかを示しています。0から1の間です。ターゲット密度に従うと、このサンプルはもっとありそうなものよりもはるかに可能性が低いので、受け入れる可能性は小さいですが、ゼロではありません。
次に進みましょう。受け入れませんでした。今度はほぼ重なっているので、おそらく受け入れるでしょう。
これをしばらく続けて、この大きなチェーンを構築しましょう。数サンプル後ではまだターゲット密度のようには見えませんが、しばらくするとこのチェーンにより、サンプルはこのターゲット密度に従い始めます。
この大きなチェーンを構築しました。チェーン内にあったすべてのサンプルをプロットすると、サンプリングしようとしていたターゲットのように見えるこの素敵なヒストグラムが得られます。驚くべきことです。
6. MCMCの実装上の課題と解決策
6.1. 初期サンプルの選択問題
Sydney Katz: 先ほど見たように、サンプル数を増やすと、望みの分布が得られます。しかし、ここに問題があります。これはすべて機能し、ターゲット分布に収束しますが、私が言った追加の条件があります。それは「無限のサンプルの極限で」ということです。残念ながら、無限のサンプルを取る時間はありません。
そのため、実際には多くのヒューリスティックやテクニックが使用されています。どのようなカーネルを使おうとも、無限のサンプルの極限では望みのものに収束します。しかし、無限のサンプルを取ろうとするのは現実的ではありません。無限ではない場合でもより良く機能するようにするためには、いくつかの工夫が必要です。
ここで、初期サンプルについて話しましょう。先ほどは失敗領域内の初期サンプルから始めました。しかし、失敗領域外の初期サンプル、例えば値が1のようなものから始めるとどうなるでしょうか。
このような初期サンプルから始めると、新しい提案サンプルもまた0の可能性を持ちます。両方とも0の可能性なので、新しいものを受け入れず、この悪い初期サンプルに固執することになります。失敗領域内の何かを提案分布から得るまで、新しいサンプルを受け入れないでしょう。
その結果、物事がうまく機能していない期間があります。分布にこのランダムな追加のピークがあり、無限サンプルの極限ではゼロサンプルのようになりますが、現時点ではそうは見えません。これは理想的ではありません。
一般に、MCMCの開始時には、開始場所によって、最初の数サンプルがサンプリングしようとしている分布のようには見えないことがあります。そのため、一般的なアプローチとして「バーンイン期間」と呼ばれるものがあります。最初の何サンプルかを捨てるだけです。例えば、最初の20サンプルを捨てるとしましょう。この領域のすべてのサンプルを捨てると、それが取り除かれます。そして残りのサンプルを取り、それ以降はすべてが合理的に見えます。
[質問：失敗分布外のものからサンプリングを続けていても、最終的には均等にならないのですか？]
質問は、失敗分布外からのサンプリングを続けても均等にならないのかということですね。はい、それが起こるのは失敗分布から始めた場合のみです。しかし一般的には、失敗に限らず、開始場所によっては最初は良いパフォーマンスが得られないことがあります。
無限のサンプルがない場合に、より良いパフォーマンスを得るためのテクニックの一つが、今話していた「バーンイン期間」です。もう一つ一般的なことは、これはマルコフ連鎖であり、次のサンプルが前のサンプルに依存しているため、サンプル間に相関があるということです。
6.2. バーンイン期間の設定
Sydney Katz: バーンイン期間は、MCMCの初期段階で生成されたサンプルを破棄する重要な技術です。先ほどお見せしたように、初期サンプルの選択によっては、最初の数回の反復で生成されるサンプルが実際のターゲット分布を代表していない可能性があります。
失敗領域外から始めた場合の例を見てみましょう。例えば値が1のサンプルから始めると、初期段階では両方のサンプル（現在のサンプルと提案されたサンプル）がターゲット密度の下で0の可能性を持ちます。つまり、新しい提案サンプルが失敗領域内に入るまで、初期サンプルに固執することになります。これは問題であり、ヒストグラムに余分なピークが現れるなど、分布の表現に悪影響を与えます。
無限のサンプル極限では、このような問題は解消されますが、私たちには無限のサンプルを取る時間がありません。そこでバーンイン期間という実用的な解決策が登場します。単純に最初のn個のサンプル（例えば最初の20個）を破棄するのです。この手法により、チェーンが実際のターゲット分布に従うようになった後のサンプルだけを保持できます。
Plutoノートブックで示したように、失敗領域外から始めると、最初はターゲット分布と大きく異なる分布が得られます。しかし、最初の一連のサンプルを破棄した後は、残りのサンプルがターゲット分布をより正確に反映するようになります。
バーンイン期間の長さを決めるのは、しばしば経験則によります。短すぎると初期サンプルの偏りが残り、長すぎると有用なサンプルが不必要に破棄されます。システムの性質やマルコフ連鎖の混合速度によっても適切なバーンイン期間の長さは変わります。
実装の観点からは、バーンイン期間を設定することは非常に簡単です。単に最初のn個のサンプルを最終結果から除外するだけです。これは後処理ステップとして行うことができ、チェーン自体の構築方法を変更する必要はありません。
バーンイン期間の設定は、サンプルの質を向上させ、結果のヒストグラムやその他の統計量がターゲット分布をより正確に反映することを保証するための重要なステップです。
6.3. サンプル間の相関と間引き
Sydney Katz: MCMCの興味深い特性の一つは、マルコフ連鎖であるため、次のサンプルが前のサンプルに依存するという点です。これは、連続するサンプル間に相関関係が生じることを意味します。この相関は、特に無限サンプルがない実際の適用では重要な課題となります。
この問題に対処するための一般的なアプローチが「間引き（thinning）」です。間引きは、チェーン内のすべてのサンプルを保持するのではなく、例えば5番目ごと、または10番目ごとのサンプルだけを保持するという手法です。これにより、保持されるサンプル間の相関を可能な限り減らすことができます。
具体的には、チェーンを構築する過程で生成されるすべてのサンプルを保存するのではなく、一定間隔でサンプルを選択します。例えば、間引き率が5の場合、5, 10, 15, ...といったインデックスのサンプルのみを最終的なサンプルセットに含めます。これは、サンプル間の依存性を減らし、より独立したサンプルのセットを得るための方法です。
間引きを行う理由は複数あります。まず、相関のあるサンプルは独立したサンプルと比較して新しい情報が少なく、効果的なサンプルサイズが小さくなります。例えば、高度に相関した1000個のサンプルは、実効的には独立した100個のサンプルと同等かもしれません。
また、相関のあるサンプルは統計的推論において問題を引き起こす可能性があります。例えば、分散の推定値が過小評価されたり、信頼区間が狭すぎたりする可能性があります。間引きを行うことで、これらの問題を軽減できます。
適切な間引き率を選択することは、システムの性質やカーネルの選択に依存します。カーネルの幅が狭い場合（つまり、新しい提案サンプルが現在のサンプルから大きく離れることがない場合）、より高い間引き率が必要になるかもしれません。逆に、カーネルの幅が広く、チェーンがより早く混合する場合、より低い間引き率で十分かもしれません。
サンプルの間引きにはトレードオフがあります。間引き率を高くすると、より独立したサンプルが得られますが、使用可能なサンプルの総数は減少します。特に計算コストが高い場合や、サンプリングに時間がかかる場合、このトレードオフは重要な考慮事項となります。
最終的には、バーンイン期間と間引きを組み合わせることで、質の高いサンプルセットを得ることができ、それによって非正規化ターゲット密度からの効果的なサンプリングが実現します。
6.4. 複数の失敗モードへの対応
Sydney Katz: MCMCを使用する際の重要な課題の一つは、複数の失敗モードが存在する場合の挙動です。このような状況をPlutoノートブックの例で見てみましょう。
ここでは、複数の失敗モードがある状況を想定しています。サンプルが低すぎる場合、そして高すぎる場合の両方で失敗とみなされます。つまり、2つの異なる失敗モードがあり、その間に失敗分布の可能性がゼロの領域があります。
このような設定でMCMCを適用した場合、何が起こるでしょうか？実際に実行してみると、一つの失敗モードからのサンプルしか得られないという問題が発生します。なぜこれが起こったのか、何が問題なのか想像できますか？
理由は、gの分布（カーネル）が現在の位置から少ししか離れない形でサンプルを提案するからです。右側の失敗モードまで移動するのが非常に難しいのです。例えば、最後の例を見てみると、現在のサンプルがあり、次のサンプルを得る可能性のある分布がここにあります。次のサンプルがこちら側（遠く離れた失敗モード）になる可能性は非常に低いです。
無限サンプルの極限では、これは最終的に起こり、この失敗モードにジャンプし、無限にサンプリングすれば戻ってくるかもしれません。しかし、再び言いますが、無限の時間がないため、より良い方法を考える必要があります。
これが「スムージング（smoothing）」と呼ばれるものの動機です。この考え方は、これらの失敗モード間を移動するためにアルゴリズムを支援する必要があるということです。それ自体ではそれを行う可能性が非常に低いからです。
失敗分布をスムーズにすることでこれを行うことができます。ゼロの位置をわずかに非ゼロにするのです。具体的な方法については、次のスムージング技術のセクションで詳細に説明します。
複数の失敗モードの問題は、特に実際のシステムで重要です。多くの実際のシステムでは、失敗は完全に異なる理由で発生する可能性があり、これらは状態空間の異なる領域に対応します。例えば、自動運転車が左に曲がりすぎて道路から外れるか、右に曲がりすぎて外れるかという失敗モードがあるかもしれません。
従来のMCMCアプローチでは、往々にして一つの失敗モードに「閉じ込められ」、他の失敗モードに十分な時間内にジャンプする可能性が極めて低くなります。これは失敗分布の不完全な表現につながり、潜在的な重大な失敗シナリオを見逃す可能性があります。
スムージングなどの技術は、このような課題に対処するために不可欠です。これにより、アルゴリズムは異なる失敗モード間をより効率的に移動でき、結果としてより包括的な失敗分布の表現が得られます。
7. スムージング技術
7.1. 失敗モード間の遷移を改善する方法
Sydney Katz: 複数の失敗モードが存在する場合、従来のMCMCでは異なる失敗モード間を移動することが困難であることを見てきました。これは、カーネルがサンプルを現在の位置から少ししか離れない場所に提案する傾向があるためです。特に、失敗モードが状態空間で遠く離れている場合、あるモードから別のモードへのジャンプは極めて低い確率でしか発生しません。
この問題を解決するためにスムージング技術を導入します。スムージングの基本的な考え方は、失敗モード間の移動をアルゴリズムが行いやすくするというものです。カーネル自体では、これらのジャンプを生成する可能性が非常に低いからです。
失敗分布をスムーズにするには、ゼロ値の位置（失敗でない領域）を少し非ゼロにします。これにより、失敗モード間の「谷」や「障壁」を越えやすくなります。
簡単に言えば、スムージングは失敗分布の「厳格さ」を緩和し、失敗領域と非失敗領域の間の遷移をより滑らかにします。これを行うことで、マルコフ連鎖はある失敗モードから次の失敗モードへとより容易に移動できるようになります。
スムージングの利点は、サンプル効率が向上し、複数の失敗モードを持つ失敗分布全体を捕捉できることです。しかし、これは純粋な失敗分布からサンプリングするという当初の目的からは少し逸脱しています。スムージングを適用すると、もはや純粋な失敗分布からサンプリングしているわけではなく、その近似からサンプリングしていることになります。
ただし、この問題にも解決策があります。スムージングした分布からサンプルを取得した後、失敗でない（つまり、元の指示関数で値が0になる）サンプルを除外することができます。これにより、スムージングの利点を活かしながらも、最終的には純粋な失敗分布からのサンプルを得ることができます。
実際には、このアプローチはスムージングした分布を提案分布として使用し、リジェクションサンプリングを適用していると考えることができます。この方法は非常に強力で、複数の失敗モードを持つ複雑なシステムにおいても効果的なサンプリングを可能にします。
次のセクションでは、このスムージング技術をより正確に定式化し、数学的に説明していきます。特に、距離関数の導入とスムージングパラメータの効果について詳細に見ていきます。
7.2. 距離関数の導入
Sydney Katz: スムージング技術をより形式的に定義するために、「距離関数」というものを導入します。この距離関数を δ(τ) と表記し、次のように定義します：τが失敗であれば0、そうでなければ0より大きい値を返す関数です。また、失敗に近づくにつれて0に近づく値を返します。
この距離関数を構築する非常に簡単な方法は、ロバスト性関数の値と0の最大値を取ることです。つまり、max(ρ(τ), 0) です。ここでρ(τ)はロバスト性関数です。成功の場合、ロバスト性は失敗に近づくにつれて低くなります。失敗の場合、ロバスト性は負の値になるので、この場合は0を返します。これは距離関数を得るための非常に簡単な方法です。
この距離関数を使って、非正規化密度を書き直すことができます。以前は、指示関数を使っていました。この指示関数は、τが仕様を満たさない（つまり失敗である）場合に1を返し、そうでなければ0を返します。ここで、失敗であるということは、この距離が0に等しいことを意味します。したがって、その条件をこの距離が0に等しいという条件に置き換えることができます。
この指示関数を考えてみると、これは実際には非常に小さな分散を持つ正規分布に似ています。距離関数δ(τ)が0のとき、指示関数は1を返します。それ以外の場合は常に0を返します。これは0を中心とした非常に小さな分散を持つ分布と似ていると言えますが、完全に同じではありません。なぜなら、今は0でない値に対しても非ゼロの値を割り当てています。
スムージングを行うために、この指示関数をこの小さな分散を持つ正規分布に置き換えることができます。こうすることで、この滑らかな分布は、失敗に対しては非ゼロの確率を割り当て、また成功に対しても非ゼロの確率を割り当てますが、失敗に近い成功ほど高い確率を持ちます。
より直感的に理解するために、これが視覚的にどのように見えるか後で示します。もう一つ注目すべき点として、このパラメータεは、どれだけスムージングを適用するかを制御します。εが0の場合、指示関数に戻ります。スムージングはありません。εが無限大に近づくと、正規分布の無限大の分散を想像できます。それは非常に広く、ほぼ一様分布のように見えます。この場合、これはp(τ)に近づきます。つまり、εが無限大の場合、名目軌道分布に近づきます。
私はこれを本当に速く説明しましたが、火曜日にも改めて説明します。今は皆さんに週末の前にこれを見せたかったのです。
多くのサンプルが失敗ではない場合がありますが、それらを拒否することで、失敗分布からのサンプルだけを残すことができます。これは、スムージングした分布を提案分布として使用し、リジェクションサンプリングを行うと考えることができます。これは非常に興味深いアプローチです。
最後に、スムージングがどれほど効果的かを示すクールなアニメーションをお見せします。これが少量のスムージングの場合です。まだ追加の失敗モードを得ていません。このスムージングされた分布が見えます。
スムージングの量を増やすと、実際に失敗モード間を行き来し始めるのが分かります。これがスムージングされた失敗分布であり、両方の失敗モードからの素晴らしいサンプルセットが得られます。
7.3. スムージングパラメータεの効果
Sydney Katz: スムージングパラメータεは、失敗分布のスムージングの程度を制御する重要なパラメータです。このパラメータの効果を理解することは、効果的なスムージング戦略を開発する上で不可欠です。
εは、距離関数δ(τ)に適用される正規分布の分散を制御します。小さいεは狭い正規分布を意味し、大きいεは広い正規分布を意味します。これがどのように働くかを見てみましょう：
εが0の場合、指示関数に戻ります。つまり、失敗の場合は1、成功の場合は0を返す関数になります。これは元の非スムージング版のMCMCに相当し、失敗モード間を移動することが難しいという問題を抱えています。
εが増加すると、失敗の周りの「確率の山」が広がります。これにより、失敗に近い成功状態にも非ゼロの確率が割り当てられます。成功状態に割り当てられる確率は、失敗状態からの距離に反比例し、距離が近いほど高い確率が割り当てられます。
εが非常に大きい（無限大に近い）場合、正規分布は非常に広くなり、ほぼ一様分布に近づきます。この極端な場合、スムージングされた分布は元の名目軌道分布p(τ)に近づきます。つまり、εが大きくなりすぎると、失敗分布の特性が薄れ、代わりに名目分布の特性が強くなります。
適切なεの値を選択することは、効果的なスムージングのために極めて重要です。値が小さすぎると、複数の失敗モード間を移動するという当初の問題が解決されません。値が大きすぎると、失敗分布の特性が失われ、名目分布に近づきすぎてしまいます。
εの最適値は、具体的な問題設定に依存します。具体的には、異なる失敗モード間の距離、失敗領域の幅、およびカーネルの性質などの要因に影響されます。実際には、いくつかの異なるε値を試し、どの値が失敗モード間の適切な移動を可能にしながらも、失敗分布の特性を維持するかを確認することが一般的です。
Plutoノートブックでの例で見たように、適切な量のスムージングを適用することで、マルコフ連鎖は複数の失敗モード間を効果的に移動できるようになります。これにより、すべての関連する失敗モードを捕捉する、より完全な失敗分布の表現が得られます。
スムージングのもう一つの利点は、失敗に近い成功状態に関する情報も提供することです。これは、システムがどのように失敗に近づくかを理解し、潜在的な脆弱性を特定するのに役立ちます。これらの「ニアミス」状態は、システムの安全性を向上させるための貴重な洞察を提供することがあります。
7.4. スムージングによる複数失敗モードのサンプリング改善
Sydney Katz: スムージング技術を適用することで、複数の失敗モードを持つシステムからのサンプリングがどのように改善されるかを具体的に見ていきましょう。Plutoノートブックのデモンストレーションは、この改善を視覚的に示しています。
まず、スムージングなしの通常のMCMCを実行すると、マルコフ連鎖は最初に遭遇した失敗モードに「閉じ込められ」てしまいます。これは、ある失敗モードから別の失敗モードへのジャンプ確率が非常に低いためです。具体的には、現在のカーネル（通常はガウス分布）が現在の位置から遠く離れた位置に新しいサンプルを提案する可能性が低いからです。
実際の例を見ると、二つの異なる失敗モードがある設定で、スムージングを適用しない場合、サンプルはすべて一方の失敗モードからのみ取得されています。これは、もう一方の失敗モードが存在するにもかかわらず、そちらにジャンプする確率が極めて低いためです。
少量のスムージング（小さなε値）を適用しても、依然として同じ問題が発生します。これは、スムージングがまだ十分でなく、失敗モード間の「谷」を埋めるのに不十分だからです。サンプルはまだ一つの失敗モードに集中しており、他のモードは捕捉できていません。
しかし、スムージング量を増やす（εを増加させる）と、状況は劇的に改善します。アニメーションで示されているように、より大きなεを使用すると、マルコフ連鎖は失敗モード間を効果的に行き来し始めます。これにより、両方の失敗モードからサンプルが取得され、より完全な失敗分布の表現が得られます。
スムージングされた分布自体も視覚化されており、失敗モード間の「谷」がどのように埋められるかを示しています。オリジナルの失敗分布では、失敗モード間にゼロ確率の領域があるのに対し、スムージングされた分布では、その領域も低いながらも非ゼロの確率を持ちます。
この改善の結果、得られるサンプルセットは両方の失敗モードをカバーしており、システムの失敗特性をより完全に表現しています。これは、安全クリティカルシステムの検証において非常に重要です。なぜなら、すべての潜在的な失敗モードを特定し、理解することが、総合的な安全対策を開発するために不可欠だからです。
スムージング技術を用いると、純粋な失敗サンプルだけでなく、失敗に近い「ニアミス」状態もサンプリングすることになります。これらのサンプルは、適切なフィルタリングプロセスを通じて除外することも、保持して追加の分析を行うこともできます。「ニアミス」状態を分析することで、システムがどのように失敗に近づくかについての貴重な洞察が得られ、より堅牢なシステム設計に役立てることができます。
このスムージング技術は、実際には私が講義の冒頭で言及した「ロバスト性」の概念に関連しています。ロバスト性関数を使用して距離関数を構築することで、失敗のバイナリな定義（失敗か成功か）から、失敗からの「距離」という連続的な尺度に移行します。これにより、失敗分布のサンプリングが大幅に改善されます。
8. 高次元システムへの応用
8.1. 1次元ガウスから振り子システムへの拡張
Sydney Katz: これまでの講義では、1次元ガウス分布という非常に単純な例を使ってリジェクションサンプリングやMCMCの概念を説明してきました。この例は理解しやすく視覚化もできますが、少し馬鹿げているようにも思えます。ガウス分布については多くのことを既に知っており、この単純なケースであれば失敗分布を解析的に計算することも可能です。これらの説明は主に直感を得るためのものでした。
しかし重要なのは、これらの手法が1次元ガウスよりも大きなシステムにも適用できるということです。例えば振り子システムの場合、50次元や100次元といった高次元の状態空間を持つシステムにスケールアップすることができます。
振り子システムへの拡張を考えると、状態は時間の経過とともに変化する位置と速度のベクトルで構成されます。各時間ステップでの状態を連結すると、非常に高次元の軌道ベクトルになります。例えば、50時間ステップのシミュレーションで、各ステップが2次元の状態（位置と速度）を持つ場合、全体の軌道は100次元になります。
このような高次元システムでは、失敗分布の視覚化は不可能ですが、サンプリング手法は依然として適用可能です。実際、これらの手法の真の価値は、解析的な計算が実行不可能な複雑なシステムにあります。
リジェクションサンプリングを振り子システムに適用する場合、提案分布の選択がより難しくなります。名目分布を使用すると、失敗が稀なシステムでは極めて非効率になる可能性があります。何千、何百万ものサンプルを生成しても、一つの失敗も見つからない可能性があるのです。
同様に、MCMCを適用する場合も、カーネルの選択が重要になります。高次元空間では、ランダムウォークの性質により、カーネルが大きすぎると新しい提案がほとんど常に拒否され、小さすぎると状態空間の探索が非常に遅くなります。これは「次元の呪い」として知られる問題です。
しかし、適切なパラメータとカーネル設計、そしてスムージング技術を用いることで、これらの手法は高次元システムでも効果的に機能します。教科書の例6.1では、振り子システムにリジェクションサンプリングを適用する際の課題と、それがどのように崩壊するかについて詳しく説明しています。
振り子システムへの拡張では、以下の点に注意する必要があります：
- 状態空間の次元が高くなるため、効率的なサンプリングがより難しくなります
- 失敗領域が高次元空間内の複雑な形状を持つ可能性があります
- 複数の失敗モードが存在し、それらが高次元空間内で遠く離れている可能性があります
- 適切なカーネルとスムージングパラメータの選択がより重要になります
それにもかかわらず、これらの手法は複雑なシステムの失敗分布をサンプリングするための強力なアプローチを提供し、安全クリティカルシステムの検証において不可欠なツールとなっています。
8.2. 高次元空間での効率的なサンプリングの課題
Sydney Katz: 高次元空間での効率的なサンプリングは、いくつかの根本的な課題を抱えています。これらの課題は「次元の呪い」として知られており、状態空間の次元が増加するにつれて顕著になります。
まず、リジェクションサンプリングを高次元システムに適用する際の主要な課題は、適切な提案分布の選択です。教科書の例6.1で示されているように、振り子システムに単純なリジェクションサンプリングを適用しようとすると、非常に低い受け入れ率に悩まされることになります。これは、高次元空間では体積が指数関数的に増加するため、失敗領域が全体の状態空間に占める割合が極めて小さくなるからです。
具体的には、名目分布を提案分布として使用した場合、50次元や100次元の振り子システムでは、失敗サンプルを1つ得るために数百万、あるいは数十億のサンプルを生成しなければならない可能性があります。これは計算的に実行不可能です。
より効率的な提案分布を手動で設計しようとしても、高次元空間では直感が働きにくく、適切な分布を見つけることが非常に困難です。また、提案分布がターゲット密度を常に上回るようにするためのスケーリング定数cの選択も、高次元では極めて難しくなります。教科書の例が示すように、保守的な値を選ぶと、cの値が非常に大きくなり、サンプリングの効率が著しく低下します。
MCMCアプローチでも、高次元空間での効率的なサンプリングには独自の課題があります。カーネルの選択が極めて重要になり、不適切なカーネルを選ぶと、連鎖が状態空間を効率的に探索できなくなります。カーネルが小さすぎると、連鎖は状態空間をゆっくりと移動し、すべての関連する領域を探索するのに非常に長い時間がかかります。逆に、カーネルが大きすぎると、ほとんどの提案サンプルが拒否され、連鎖が「立ち往生」する可能性があります。
複数の失敗モードが存在する場合、高次元空間ではこれらのモード間の距離がさらに大きくなる傾向があり、標準的なMCMCでは一方のモードから他方へ移動することが実質的に不可能になります。この問題に対処するためには、より洗練されたカーネル設計やスムージング技術が不可欠です。
スムージング技術は高次元システムでも効果的ですが、適切なεの値の選択はさらに難しくなります。εが小さすぎると失敗モード間の移動が不十分になり、大きすぎると失敗分布の特性が失われます。
これらの課題に対処するためのアプローチには以下が含まれます：
- 適応型MCMC法：連鎖の進行に基づいてカーネルを動的に調整する方法
- ハミルトニアンモンテカルロ法：勾配情報を利用して、より効率的に状態空間を探索する方法
- 並列テンパリング：複数の温度パラメータを持つ複数の連鎖を並行して実行し、連鎖間でサンプルを交換する方法
- 次元削減技術：問題の本質的な次元を減らすことで、サンプリングを効率化する方法
高次元システムでの効率的なサンプリングは依然として活発な研究分野であり、安全クリティカルシステムの検証にこれらの技術を適用する際には、問題の特性に応じた慎重なアプローチが求められます。しかし、適切に設計された場合、これらの手法は複雑なシステムの失敗分布を理解する上で非常に価値のある情報を提供することができます。
