【確率論】確率分布まとめ

Poisson分布

定義

 p(k; \lambda)= e^{-\lambda}\frac{\lambda^k}{k!}

導出

Poisson分布は二項分布の np=\lambda, n \rightarrow \infty(p \rightarrow 0) にして定義されている.
二項分布

\begin{aligned}
b(k; n,p) 
&= \dbinom{n}{k}p^k(1-p)^{n-k}\\
&= \frac{n!}{k!(n-k)!}p^k(1-p)^{n-k}\\
&= \frac{n(n-1)\cdots (n-k+1)}{k!}p^k(1-p)^{n-k}\\
&= \frac{n(n-1)\cdots (n-k+1)}{n \cdot n \cdots n}n^kp^k\frac{1}{k!}(1-p)^{n-k}\\
&= \frac{n(n-1)\cdots (n-k+1)}{n \cdot n \cdots n}(np)^k\frac{1}{k!}\frac{(1-p)^{n}}{(1-p)^{k}}\\
&= \frac{n(n-1)\cdots (n-k+1)}{n \cdot n \cdots n}\lambda^k\frac{1}{k!}\frac{((1-\frac{\lambda}{n})^{-\frac{n}{\lambda}})^{-\lambda}}{(1-p)^{k}}\\
& \rightarrow 1\cdot \lambda^k \cdot \frac{1}{k!} \cdot e^{-\lambda}\\
&= e^{-\lambda}\frac{\lambda^k}{k!}
\end{aligned}

積分

 f(x)=e^xマクローリン展開 f^{(n)}(x)=e^xより f^{(n)}(0)=1であるから,
 e^x=f(0)+f'(0)+\frac{1}{2!}f''(0)x^2+\cdots + \frac{1}{n!}f^{(n)}x^n+\cdots =1+x+\frac{1}{2!}x^2+\cdots + \frac{1}{n!}x^n+\cdots

\begin{aligned}
\sum^{\infty}_{k=0}e^{-\lambda}\frac{\lambda^k}{k!}
&= e^{-\lambda}\sum_{k=0}^{\infty}\frac{\lambda^k}{k!}\\
&= e^{-\lambda}\left(\frac{\lambda^0}{0!}+\frac{\lambda^1}{1!}+\frac{\lambda^2}{2!}+\cdots + \frac{\lambda^n}{n!}+\cdots \right)\\
&= e^{-\lambda}e^{\lambda}\\
&= 1
\end{aligned}

期待値


\begin{aligned}
E(X)
&= \sum_{k=0}^{\infty}ke^{-\lambda}\frac{\lambda^{k}}{k!}\\
&= e^{-\lambda}\sum_{k=0}^{\infty}k\frac{\lambda^{k}}{k!}\\
&= e^{-\lambda}\sum_{k=1}^{\infty}k\frac{\lambda^{k}}{k!}\\
&= e^{-\lambda}\lambda\sum_{k=1}^{\infty}\frac{\lambda^{k-1}}{(k-1)!}\\
&= e^{-\lambda}\lambda\sum_{m=0}^{\infty}\frac{\lambda^{m}}{m!}\\
&= e^{-\lambda}\lambda e^{\lambda}\\
&= \lambda
\end{aligned}

分散

 k^2=k(k-1)+kを利用する.

\begin{aligned}
E(X^2)
&=\sum_{k=0}^{\infty}k^2e^{-\lambda}\frac{\lambda^k}{k!}\\
&= \left(\sum_{k=0}^{\infty}k(k-1)\frac{\lambda^{k-2}}{k(k-1)(k-2)!}\lambda^2+e^{\lambda}E(X)\right)\\
&= e^{-\lambda}\left(\sum_{k=2}^{\infty}\frac{\lambda^{k-2}}{(k-2)!}\lambda^2+e^{\lambda}\lambda\right)\\
&= e^{-\lambda}\left(\sum_{m=0}^{\infty}\frac{\lambda^{m}}{m!}\lambda^2+e^{\lambda}\lambda\right)\\
&= e^{-\lambda}(e^{\lambda}\lambda^2+e^{\lambda}\lambda)\\
&= \lambda^2 + \lambda\\
V(X) 
&= E(X^2) -E(X)^2\\
&= \lambda^2 + \lambda - \lambda^2\\
&= \lambda
\end{aligned}

import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as stats
from IPython.core.pylabtools import figsize
figsize(12.5, 4)

poi = stats.poisson # Poisson分布
lambda_ = [1.5, 5.5, 12] # Poisson分布のパラメタ(強度λ)
colors = ["#348ABD", "#A60628", "#66FFB2"]

a = np.arange(16) # [0, 1, 2, ..., 15]
# poi.pmf(a, λ)でPoisson分布の値を得る
plt.bar(a, poi.pmf(a, lambda_[0]), color=colors[0], label="$\lambda = %.1f$"%lambda_[0], alpha=0.6, edgecolor=colors[0], lw=1)
plt.bar(a, poi.pmf(a, lambda_[1]), color=colors[1], label="$\lambda = %.1f$"%lambda_[1], alpha=0.6, edgecolor=colors[1], lw=1)
plt.bar(a, poi.pmf(a, lambda_[2]), color=colors[2], label="$\lambda = %.1f$"%lambda_[2], alpha=0.6, edgecolor=colors[2], lw=1)
plt.legend()
plt.ylabel("Probability of $k$") #kの確率
plt.xlabel("$k$")

f:id:umashika5555:20200228004910p:plain

幾何分布

定義

 P\left\{X=k\right\}=(1-p)^{k}p

意味

ベルヌーイ試行で初めて成功するまでの試行回数の確率

期待値

 
\begin{aligned}
E(X)
&= \sum_{k=0}^{\infty}k(1-p)^kp\\
&=p\lim_{n\rightarrow \infty}\sum_{k=0}^{n}k(1-p)^{k}\\
& =p\lim_{n\rightarrow \infty} \frac{1-p}{p^2}\left(1-(1-p)^n(1+np^2)\right)\\
&= p \frac{1-p}{p^2}\\
&= \frac{1-p}{p}
\end{aligned}

分散

 
\begin{aligned}
E(X^2)
&=\sum^{\infty}_{k=0}k^2(1-p)^{k}p\\
&=p\sum^{\infty}_{k=0}k^2(1-p)^k\\
&=p\lim_{n\rightarrow \infty}\sum^{n}_{k=0}k^2(1-p)^k\\\\
&= p \lim_{n\rightarrow \infty}\left(\frac{1}{p}\left\{2\sum^{n}_{k=1}k(1-p)^k - \sum_{k=1}^n(1-p)^k\right\}\right)\\
&= p \lim_{n\rightarrow \infty} \frac{1}{p}\left[\left(\frac{2}{p}-1\right)\frac{(1-p)\{1-(1-p)^n\}}{p}-\frac{2n}{p}(1-p)^n\right]\\
&= p\frac{1}{p}\left(\frac{2}{p}-1\right)\frac{1-p}{p}\\
&= \frac{p^2-3p+2}{p^2}\\
V(X)
&= E(X^2) - E(X)^2\\
&= \frac{p^2-3p+2}{p^2}-\left(\frac{1-p}{p}\right)^2\\
&= \frac{1-p}{p^2}
\end{aligned}

scipy.stats.geomでは (1-p)^{k-1}pで定義されているため, k=0は定義なし.
成功確率が高いときは大きい試行回数の値で確率が小さくなり, 反対に成功確率が小さいときは大きい試行回数でも確率が比較的に小さくならない.

import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as stats
from IPython.core.pylabtools import figsize
figsize(12.5, 4)

geo = stats.geom # Geometric分布
p = [0.1, 0.5, 0.9] # Geometric分布のパラメタ(確率p)
colors = ["#348ABD", "#A60628", "#66FFB2"]

k = np.arange(16) # [0, 1, 2, ..., 15]
# geo.pmf(k,p)でGeometric分布の値を得る
plt.bar(k, geo.pmf(k, p[0]), color=colors[0], label="$p = %.1f$"%p[0], alpha=0.6, edgecolor=colors[0], lw=1)
plt.bar(k, geo.pmf(k, p[1]), color=colors[1], label="$p = %.1f$"%p[1], alpha=0.6, edgecolor=colors[1], lw=1)
plt.bar(k, geo.pmf(k, p[2]), color=colors[2], label="$p = %.1f$"%p[2], alpha=0.6, edgecolor=colors[2], lw=1)

plt.xticks(k, k)
plt.legend()
plt.ylabel("Probability of $k$") #kの確率
plt.xlabel("$k$")

f:id:umashika5555:20200228010632p:plain

Gauss分布

定義

 N(x ; \mu, \sigma^2)=\frac{1}{\sqrt{2\pi\sigma^2}}\exp{\left(-\frac{(x-\mu)^2}{2\sigma^2}\right)}

積分

まずガウス積分 \int_{-\infty}^{+\infty}\exp{(-x^2)}dx=\sqrt{\pi}の導出を行う.
 x=r \cos \theta , y=r \sin \thetaと置き, ヤコビアン rと計算できる.

\begin{aligned}
I =  \int_{-\infty}^{+\infty}\exp{(-x^2)}dx 
\end{aligned}
とおく.

\begin{aligned}
I^2 
&=  \int_{-\infty}^{+\infty}\exp{(-x^2)}dx \int_{-\infty}^{+\infty}\exp{(-x^2)}dx \\
&=  \int_{-\infty}^{+\infty}\int_{-\infty}^{+\infty} \exp{(-(x^2+y^2))}dxdy \\
&= \int_{0}^{2\pi}d\theta \int^{\infty}_{0} \exp{(-r^2)}r dr\\
&= 2\pi(-\frac{1}{2})\left[\exp{(-r^2)}\right]^{\infty}_{0}\\
&= 2\pi(-\frac{1}{2})(0-1)\\
&=\pi\\
\therefore I &=\sqrt{\pi}\ (\because\ I>0)
\end{aligned}
次に,  \int_{-\infty}^{+\infty}\exp{(-ax^2)}dx=\sqrt{\frac{\pi}{a}}の導出を行う.
これは \sqrt{a}x = tと置きガウス積分に代入できる.
その際,  dx = \frac{1}{\sqrt{a}}dtよりガウス積分の結果 \sqrt{\pi}に係数 \frac{1}{\sqrt{a}}が付き,
最終的な結果が \sqrt{\frac{\pi}{a}}となる.

これらの結果から, 正規分布積分が計算できる.
 
\begin{aligned}
\int^{+\infty}_{-\infty}\frac{1}{\sqrt{2\pi\sigma^2}}\exp{\left(-\frac{(x-\mu)^2}{2\sigma^2}\right)}dx
&= \frac{1}{\sqrt{2\pi\sigma^2}} \int^{+\infty}_{-\infty}\exp{\left(-\frac{(x-\mu)^2}{2\sigma^2}\right)}dx\\
&= \frac{1}{\sqrt{2\pi\sigma^2}} \sqrt{\frac{\pi}{\frac{1}{2\sigma^2}}}\\
&= \frac{1}{\sqrt{2\pi\sigma^2}} \sqrt{2\pi\sigma^2}\\
&= 1
\end{aligned}

期待値

 x=(x-\mu) + \muとして, 後ろの項は定数なので, 確率分布のみの積分となり, 確率分布の性質から積分して1になる.
 
\begin{aligned}
E(X) 
&= \int_{-\infty}^{\infty}xf(x)dx\\
&= \int_{-\infty}^{\infty}(x-\mu)f(x)dx + \mu\int_{-\infty}^{\infty}f(x)dx\\
&= 0 + \mu \cdot 1\\
&= \mu
\end{aligned}

分散

 x^2 = (x-\mu)^2 + 2\mu x - \mu^2を用いる.

\begin{aligned}
E(X^2)
&= \int_{-\infty}^{\infty}x^2N(x; \mu, \sigma^2)dx\\
&= \frac{1}{\sqrt{2\pi\sigma^2}}\int_{-\infty}^{\infty}x^2\exp\left\{-\frac{(x-\mu)^2}{2\sigma^2}\right\}dx\\
&= \frac{1}{\sqrt{2\pi\sigma^2}}\int_{-\infty}^{\infty}(x-\mu)^2\exp\left\{-\frac{(x-\mu)^2}{2\sigma^2}\right\}dx + 2\mu \int_{-\infty}^{\infty}x\frac{1}{\sqrt{2\pi\sigma^2}}\exp\left\{-\frac{(x-\mu)^2}{2\sigma^2}\right\}dx - \frac{\mu^2 }{\sqrt{2\pi\sigma^2}}\int_{-\infty}^{\infty}\exp\left\{-\frac{(x-\mu)^2}{2\sigma^2}\right\}dx\\
&= \frac{1}{\sqrt{2\pi\sigma^2}}\int_{-\infty}^{\infty}(x-\mu)^2\exp\left\{-\frac{(x-\mu)^2}{2\sigma^2}\right\}dx + 2\mu E(X) - \frac{\mu^2 }{\sqrt{2\pi\sigma^2}}\int_{-\infty}^{\infty}\exp\left\{-\frac{(x-\mu)^2}{2\sigma^2}\right\}dx\\
\end{aligned}
第1項は,  (x-\mu)=tとして,  e^{-t^2} = -\frac{(e^{-t^2})'}{2t}と見なし部分積分して \sigma^2
第2項は,  2\mu E(X) = 2\mu\cdot\mu = 2\mu^2\
第3項はガウス積分を行い,  -\frac{\mu^2}{\sqrt{2\pi\sigma^2}}\cdot\sqrt{\frac{\pi}{\frac{1}{2\sigma^2}}}= -\frac{\mu^2}{\sqrt{2\pi\sigma^2}}\cdot \sqrt{2\pi\sigma^2}=-\mu^2
よって,  E(X^2) = \sigma^2 + 2\mu^2 - \mu^2 = \sigma^2 +\mu^2
したがって分散は,  V(X) = E(X^2) - E(X)^2 = \sigma^2 + \mu^2 -\mu^2 = \sigma^2

積率母関数


\begin{aligned}
\int_{-\infty}^{\infty}e^{\theta x}N(x; \mu, \sigma^2)dx
&= \frac{1}{\sqrt{2\pi\sigma^2}}\int_{-\infty}^{\infty}\exp{\left(\theta x-\frac{(x-\mu)^2}{2\sigma^2}\right)}dx
\end{aligned}
ここで exp()内をxについて整理する.

\begin{aligned}
g(x) &:= \theta x - \frac{(x-\mu)^2}{2\sigma^2}\\
&= -\frac{1}{2\sigma^2}\left(x-(\sigma^2\theta+\mu)\right)^2+\frac{\sigma^2\theta^2}{2}+\mu\theta
\end{aligned}
この形はガウス積分なので係数と合わせて1となる.

\begin{aligned}
e^{\frac{\sigma^2\theta^2}{2}+\mu\theta} \cdot \frac{1}{\sqrt{2\pi\sigma^2}}\int_{-\infty}^{\infty}\exp{\left(-\frac{\left(x-(\sigma^2\theta+\mu)\right)^2}{2\sigma^2}\right)}dx
&= e^{\frac{\sigma^2\theta^2}{2}+\mu\theta} \cdot 1\\
&= e^{\frac{\sigma^2\theta^2}{2}+\mu\theta}
\end{aligned}

多次元Gauss分布

 

\begin{aligned}
N(\boldsymbol{x}; \boldsymbol{\mu}, \boldsymbol{\Sigma}) = \frac{1}{(2\pi)^{\frac{D}{2}}}\frac{1}{|\Sigma|^{\frac{1}{2}}}\exp\left\{-\frac{1}{2}(\boldsymbol{x}-\boldsymbol{\mu})^{T}\boldsymbol{\Sigma}^{-1}(\boldsymbol{x}-\boldsymbol{\mu})\right\}
\end{aligned}

指数分布

定義

 f(x) = 
  \begin{cases}
    \lambda e^{-\lambda x}\ &(x\geq 0)\\
    0\ &(x<0)
  \end{cases}

期待値


\begin{aligned}
E(X)
&= \int_{-\infty}^{\infty}xf(x)dx\\
&= \int_{0}^{\infty}x\lambda e^{-\lambda x}dx\\
&= \lambda \left(-\frac{1}{\lambda}\right)\int_{0}^{\infty}x(e^{-\lambda x})'dx\\
&= - \int_{0}^{\infty}x(e^{-\lambda x})'dx\\
&= -[xe^{-\lambda x}]^{\infty}_{0} + \int_{0}^{\infty}e^{-\lambda x}dx\\
&= (0-0) + [-\frac{1}{\lambda}e^{-\lambda x}]^{\infty}_{0}\\
&= \frac{1}{\lambda}
\end{aligned}

分散

 
\begin{aligned}
E(X^2)
&= \int_{-\infty}^{\infty}x^2f(x)dx=\lambda\int^{\infty}_{0}x^2e^{-\lambda x}dx\\
&= \lambda \left(\frac{1}{\lambda}\right)\int^{\infty}_{0}x^{2}(e^{-\lambda x})'dx\\
&=-\int^{\infty}_{0}x^{2}(e^{-\lambda x})'dx\\
&= -[x^2(e^{-\lambda x})]^{\infty}_{0} + 2\int^{\infty}_{0}xe^{-\lambda x}dx\\
&= 0 + 2\left(\frac{E(X)}{\lambda}\right)\\
&= \frac{2}{\lambda^2}\\
V(X)
&= E(X^2)-E(X)^2\\
&= \frac{2}{\lambda^2} - \left(\frac{1}{\lambda}\right)^2\\
&= \frac{1}{\lambda^2}
\end{aligned}

import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as stats
from IPython.core.pylabtools import figsize
figsize(12.5, 4)

expo = stats.expon
lambda_ = [0.5, 1, 1.5] # 指数分布のパラメタλ
colors = ["#348ABD", "#A60628", "#66FFB2"]

a = np.linspace(0, 4, 100) # 本来は連続値だがグラフにする上で離散値をつくる

for l, c in zip(lambda_, colors):
  # 指数分布は scale=1/λ を入力する
  plt.plot(a, expo.pdf(a, scale=1. / l), lw=3, color=c, label="$\lambda = %.2f$"%l)
  plt.fill_between(a, expo.pdf(a, scale=1. /l), color=c, alpha=.33)
  

plt.legend()
plt.ylabel("$f(k$)") #kの確率
plt.xlabel("$k$")

f:id:umashika5555:20200228005623p:plain

二項分布

定義

 b(k; n, p) = \dbinom{n}{k}p^k(1-p)^{n-k}
k: ベルヌーイ試行についてある事象が成功する回数
n: すべての試行回数
p: ベルヌーイ試行についてある事象が成功する確率

積分

二項定理より自明
 \sum_{k=0}^{n}b(k: n,p) = \sum_{k=0}^{n}\dbinom{n}{k}p^k(1-p)^{n-k}=\left(p+(1-p)\right)^n=1^n=1

期待値


\begin{aligned}
E(X)
&=\sum_{k=0}^{n}kb(k;n,p)\\
&=0+\sum_{k=1}^{n}k\dbinom{n}{k}p^{k}(1-p)^{n-k}\\
&=\sum_{k=1}^{n}\frac{n(n-1)\cdots (n-k+1)}{(k-1)!}p^k(1-p)^{n-k}\\
&= np \sum_{k=1}^{n}\frac{(n-1)\cdots (n-k+1)}{(k-1)!}p^{k-1}(1-p)^{n-k}\\
&= np \sum_{l=0}^{n-1}\frac{(n-1)\cdots (n-l)}{l!}p^{l}(1-p)^{(n-1)-l}\\
&= np(p+(1-p))^{n-1}\\
&= np
\end{aligned}

分散

 V(X) = E(X^2) - E(X)^2を用いる.
 E(X^2)の計算では, 二項定理の係数を消すためと, 期待値の定義式に落とし込むため,
 k^2 = k(k-1) +k という処理を施す. 前半は, 期待値の計算のように二項定理の係数項が2つ分うち消える.
後半は, 期待値の定義の式と一致する.

\begin{aligned}
E(X^2)
&=\sum_{k=0}^{n}k^2b(k;n,p)\\
&=\sum_{k=2}^{n}k(k-1)\dbinom{n}{k}p^{k}(1-p)^{n-k} + \sum_{k=1}^{n}k\dbinom{n}{k}p^{k}(1-p)^{n-k}\\
&=n(n-1)\sum_{k=2}^{n}\frac{(n-2)\cdots (n-k+1)}{(k-2)!}p^k(1-p)^{n-k}+np\\
&=n(n-1)\sum_{m=0}^{n-2}\frac{(n-2)\cdots (n-m-1)}{m!}p^{m+2}(1-p)^{n-m-2}+np\\ 
&=p^{2}n(n-1)\sum_{m=0}^{n-2}\frac{(n-2)\cdots ((n-2)-m+1)}{m!}p^{m}(1-p)^{(n-2)-m}+np\\
&=p^{2}n(n-1)\left(p+(1-p)\right)^{n-2}+np\\
&=p^{2}n(n-1)+np\\
&=n^2p^2-np(1-p)\\

V(X)
&=E(X^2)-E(X)^2\\
&=n^2p^2+np(1-p)-(np)^2\\
&=np(1-p)\\
\end{aligned}

n=10のとき, p=0.1だと1回成功する確率が38%くらい. p=0.5だと5回成功する確率が25%くらい., p=0.9だと9回成功する確率が38%くらいと読み取れる.
実際, 代入して計算すると,
 b(k=1; n=10, p=0.1)=0.3874204890000005
 b(k=5, n=10, p=0.5)=0.24609375000000025
 b(k=9; n=10, p=0.9)=0.38742048900000037

import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as stats
from IPython.core.pylabtools import figsize
figsize(12.5, 4)

binom = stats.binom # Binomial分布
n_ = [5, 10]
p_ = [0.1, 0.5, 0.9] # Binomial分布のパラメタ(確率p)
colors = ["#348ABD", "#A60628", "#66FFB2"]

for n in n_:
    for i, p in enumerate(p_):
        k = np.arange(n+1) # [0, 1, 2, ..., 15]
        # binom.pmf(k, n ,p)でBinomial分布の値を得る
        plt.bar(k, binom.pmf(k, n, p), color=colors[i], label="$n= %d, p = %.1f$"%(n,p), alpha=0.6, edgecolor=colors[i], lw=1)
  
    plt.xticks(k, k)
    plt.legend()
    plt.ylabel("Probability of $k$") #kの確率
    plt.xlabel("$k$")
    plt.show()
    plt.close()

f:id:umashika5555:20200228012208p:plain
f:id:umashika5555:20200228012220p:plain