機械と学習する

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

トップ > 確率分布 > 逆関数法を利用して切断指数分布に従う乱数を取得する

逆関数法を利用して切断指数分布に従う乱数を取得する

確率分布 統計モデリング

【概要】

  • ちょっと切断指数分布に従う乱数が必要になったので、逆関数法を使って乱数を生成する方法をまとめた

【目次】

はじめに

切断指数分布に従う乱数の生成が必要になったのですが、ちょっと調べたところ、切断指数分布に従う乱数を生成するためのAPIはあまり用意されてない？っぽかったです。なので本記事では、逆関数法を使って切断指数分布に従う乱数の生成についてまとめます。

（この記事を書く中でさらに調べたところscipyにtruncexponというAPIが用意されていました。なので、切断指数分布が必要な方はこちらを利用したら良いと思います。。。）

本記事に興味がある方はこちらも参考にしてください。

learning-with-machine.hatenablog.com

【トップに戻る】

切断分布

切断分布とは、値域が切断された確率密度関数です。

正規分布を例にするとこんな感じ。 オレンジ線が\mathcal{N}(x | \mu=0.0, \sigma^2=1.0)の確率密度関数。青線が同じパラメータで-1から1までの領域に限定した切断正規分布。

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

領域が切断されているため、積分して1になるという密度関数の定義を満たすために、値の定義される領域で密度関数の値は切断分布の方が大きくなっています。なので、単純に密度関数を切るだけではないということに注意が必要です。

切断正規分布などは稀に使われているところを見かけます。例えば、事前知識として値の領域が限定されていることがわかっているパラメータの推論をする際の事前分布とか。

【トップに戻る】

逆関数法

逆関数法とは、累積密度関数をの逆関数を利用して、一様分布に従う乱数（一様乱数）から所望の確率密度関数に従う乱数を生成する方法です。

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

上図を見れば一目瞭然で、乱数を得たい確率密度関数をf(x)として、その累積分布関数F(x)を考えます。F(x)は0から1の範囲の関数なので(確率の定義から)、一様乱数に従うUをF(x)の逆関数F^{-1}(U)に入力することで、f(x)に従う乱数xに変換することができます。


\begin{align}
 X = F^{-1}(U)
\end{align}

ということで、逆関数法を利用するには累積分布関数の逆関数F^{-1}(U)を導出する必要がありますが、複雑な関数の場合は解析的に導出することが難しい場合もあります。 この点については、僕が以前書いたブログで、ノンパラメトリックな分布に対して逆関数法を適用してサンプルを取得している例があります。

【トップに戻る】

指数分布と切断指数分布

ここまでで逆関数法がわかったので、指数分布と切断指数分布にそれぞれ従う乱数を取得してみます。

指数分布

指数分布は以下の式で定義されています。


\begin{align}
 f(x) = \mu \exp(-\mu x)
\end{align}

分布関数は、これを積分すればよく、以下の形として知られています。


\begin{align}
 F(x) = 1 - \exp(-\mu x)
\end{align}

逆関数F^{-1}(\cdot)は、上記のF(x)をxについて解きます。


\begin{align}
 e^{-\mu x} &= 1 - U \\\
 x &= - \frac{1}{\mu} \log(1-U)
\end{align}

これで一様分布に従うUから指数分布\mathrm{Exp}(x | \mu)に従う乱数を取得できます。

切断指数分布

同様にして切断指数分布を考えます。 切断点をTとして、xが0からTまでの範囲は\exp(-\mu x)に比例します。それ以外の領域では0となります。 なので、0からTの範囲で積分して正規化定数を導出します。


\begin{align}
 \int^T_0 e^{-\mu x}dx = \frac{1-e^{-\mu T}}{\mu} \\\
 f(x) = \frac{\mu}{1-e^{-\mu T}} e^{-\mu x}
\end{align}

累積分布関数は


\begin{align}
 F(x) &= \int^x_0 \frac{\mu}{1-e^{-\mu T}} e^{-\mu x} dx \\\
  &= \frac{1}{1-e^{-\mu T}} \left( 1-e^{-\mu x} \right)
\end{align}

これでF(x)をUとして、xについて解くと


\begin{align}
 1-e^{-\mu x} &= (1-e^{-\mu T}) U \\\
 -\mu x &= \log\left\{ 1- (1-e^{-\mu T}) U \right\} \\\
 x &= -\frac{1}{\mu} \log\left\{ 1- (1-e^{-\mu T}) U \right\}
\end{align}

ということで、切断指数分布に従う乱数を得るための変換式を導出できました。

【トップに戻る】

試してみる

詳細は以下のnotebookを参照してください。

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

こんな感じで、それぞれの分布に従う乱数を生成できました。

「はじめに」でも書きましたが、scipyにはtrancexponという関数が用意されているみたいです。ということで、scipy使えば逆関数法による導出は要らなかった。。。