機械と学習する

統計解析、機械学習について学習したことをまとめていきます

トップ > 確率分布 > 数値的逆関数法でノンパラメトリックな確率分布からサンプルを取得する

数値的逆関数法でノンパラメトリックな確率分布からサンプルを取得する

確率分布 統計モデリング

【概要】

  • ちょっと必要になったので、ノンパラメトリックに推定した確率分布からのサンプル取得を実装した

【目次】

はじめに

現実のデータ（赤外線距離センサのデータなど）は、環境要因などでノイズが大きかったり、特殊な構造のノイズが乗っていたりすることが多いと思います。 そのような複雑な構造のデータを分析するにあたって、混合モデルなどを駆使すればだいたい既知の確率分布を利用してモデルを構築することはできます。

しかし、現実的にあまり複雑な構造のデータを頑張ってモデリングするのはきついときがあります。モデル構造を頑張って考えるよりもまずは動かしてみたいという要求も現実には多いと思います。

そこで、既知の確率分布を仮定しないノンパラメトリックな確率密度推定法（KDEなど）を利用することで、確率密度関数の推定が行えます。

ただ、密度関数の推定ができただけでは、確率モデルとして組み込んで利用することが難しいです。推定した密度関数からサンプルの取得などが行える必要があります。

ということで本記事では、ノンパラメトリックに推定した密度関数から数値的逆関数法を利用してサンプルの取得をやってみようと思います。

【トップに戻る】

ノンパラメトリックな確率密度推定

確率密度の推定はざっくりと分類して、既知の確率分布にデータを当てはめる「パラメトリックな手法」とデータだけから数値的に確率密度を推定する「ノンパラメトリックな手法」の二つがあります（参考文献[1]）。

前者はデータ構造に正規分布を仮定して平均や分散のパラメータを推論する方法が代表的です。また、ガウス混合モデルはガウス分布を多数混合させることで任意の密度関数を任意の精度で近似できるとされているようです（参考文献[1]）。

後者はヒストグラムを計算するイメージです。KDE(カーネル密度推定)などが代表的な手法です（参考文献[1]）。KDEとヒストグラム密度推定(HDE)の実装は当ブログでも扱いました。

learning-with-machine.hatenablog.com

今回はノンパラメトリックに推定した密度関数からサンプルを取得するということを扱います。

【トップに戻る】

数値的逆関数法によるサンプルの取得

KDEなどで密度関数を推定すると、任意の点での確率密度の推定値を得ることができます。しかしこれだけでは、確率モデルとして扱いにくいです。MCMCなど数値的な推論を行うには確率モデルからサンプルを取得することが必要になります。

一般に、正規分布に従う乱数の生成なども一様乱数と逆関数法を利用して生成されます。

逆関数法

既知の確率分布p(x)に従う乱数を生成することを考えます。

p(x)の累積分布関数F(v)


\begin{align}
 F(v) = \int^v_{-\infty} p(x)dx
\end{align}

この累積分布関数の逆関数v=F^{-1}(u)を導くことができれば、一様乱数uからp(x)に従う乱数vが得られます。

数値的逆関数法

ということなんですが、ここで問題になるのは、F^{-1}(\cdot)を解析的に求められないということです（p(x)が数値的に与えられている）。

そこで、F^{-1}(\cdot)を離散化して近似します。手順は次の通りです（参考文献[2]など）。

  1. データから確率分布p(x)をなんとかして推定する（ヒストグラム近似やKDEなど）
  2. p(x)の領域を適当にJ区間に分割（幅は一定でなくてもOK）
  3. 離散化したp(x)に基づいて累積分布関数F(x)を算出
  4. 分布関数の逆関数F^{-1}(\cdot)を定義

【確率分布p(x)の推定】

ノンパラメトリックに推定する場合は、KDEなどが代表的です。

【p(x)の領域を離散化】

確率分布p(x)の領域をx_{min} \sim x_{max}の区間と定義し、その中をJ区間に分割します。それぞれの分割位置をv_j、分割領域の幅をdv_j = v_j - v_{j-1}とします。

【累積分布関数F(x)の算出】

分布関数を分割区間毎に累積して算出します。


\begin{align}
 F_j = \frac{\sum^{j}_{l=1}p(v_l)dv_l }{\sum^{J}_{l=1}p(v_l)dv_l}
\end{align}

分母は正規化のための定数ですね。

【逆関数F^{-1}(\cdot)を定義】

F_jとして離散化した分布関数が定義できたので、これの逆関数を次のように定義します。

  • [0, 1]の範囲のuが、[tex:F{l-1} < u < F{l}]となるlを探す
  • v = \frac{v_{l-1} + v_{l}}{2}とする

この操作のイメージを図にしてみました。

f:id:hippy-hikky:20210218105728p:plain

青線が累積密度(F(v))を導出する流れですが、uからvを求めたいので向きが逆になります。

これで一様乱数uから任意の確率分布p(x)に従う乱数をサンプルすることができます。

【トップに戻る】

実装

実装して確認してみます。

トイデータとしてベルヌーイ分布から適当に乱数を生成しました。ベルヌーイ分布を選定した理由は定義域が[0, 1]で範囲が決まっているからという理由だけです。

f:id:hippy-hikky:20210218111954p:plain
サンプルデータ. 実線はKDEでの推定密度関数．

このデータからKDEで推定した密度関数を使ってサンプルを取得してみます。

f:id:hippy-hikky:20210218112108p:plain
数値的逆関数法でのサンプル結果．サンプル点のヒストグラム，KDEでの推定密度関数(オレンジ線), 設定した密度関数(青線)．

想定通りにサンプルが取得できていそうなことがわかりますね。

次に、もうちょっと複雑な例として、ベルヌーイ分布の混合分布を試します。サンプルデータは次の通りです。

f:id:hippy-hikky:20210218112356p:plain
サンプルデータその２．

こちらも同様に推定した密度関数に基づいてサンプルしてみましょう。

f:id:hippy-hikky:20210218112449p:plain
数値的逆関数法でのサンプル結果．

こちらも想定通りですね。混合分布ですが、明示的に混合分布を定義するのではなく、データだけで混合分布からのサンプルが行えていることがわかります。

実験コードは次節で添付しているので興味あるかたは参考にしてもらえたらと思います。効率化は考えていないので、計算は重いです。

実装コード

数値的に分布関数の導出とサンプルを実行するEmpiricalDistクラスは次のように実装しました。

【トップに戻る】

おわりに

ということで今回は、ノンパラメトリックな確率密度推定からサンプルの取得までを実装して確認しました。

ノンパラメトリックで確率分布を推定するので、混合分布とか、データ構造が複雑なデータもあまり深く考えずに利用できそうです。既知の確率分布が当てはめられるようにすることで計算が速くできるはずなので、使いどころは考える必要がありますけどね。

あと、ソースコード添付してますが、効率化を考えないと重くてなかなか辛いです。

【トップに戻る】

参考資料

【トップに戻る】