コーシー分布とラプラス分布

コーシー分布

$${\cal{N}(0,1)}$$の独立な確率変数$${z,\ z'}$$の比である$${x=\displaystyle{\frac{z}{z'}}}$$が従う確率分布を、最頻値$${0}$$、半値巾$${1}$$の標準コーシー分布と呼び、確率密度関数は$${f(x)=\displaystyle{\frac{1}{\pi(x^2+1)}}}$$で与えられる。
これを一般化して、最頻値$${a}$$、半値巾$${b}$$としたコーシー分布の確率密度関数は、
$${f(x)=\displaystyle{\frac{1}{\pi}\cdot\frac{b}{(x-a)^2+b^2}}}$$
と与えられる。
期待値は、
$${E[x]=\displaystyle{\frac{1}{\pi}\int^\infty_{-\infty}\frac{bx}{(x-a)^2+b^2}=\frac{1}{\pi}\int^\infty_{-\infty}dx\Big( \frac{b(x-a)}{(x-a)^2+b^2}+\frac{ba}{(x-a)^2+b^2}\Big) }}$$
　第二項において、$${x+a=b\tan x}$$と置換積分すれば
$${E[x]=\displaystyle{ \frac{1}{\pi}\Big\{ \Big[\frac{b}{2}\log((x-a)^2+b^2)\Big]^\infty_{-\infty} +a[\tan^{-1}x]^\infty_{-\infty}\Big\} }}$$
$${\displaystyle{=\frac{1}{\pi}\Big\{ \lim_{\alpha\to\infty,\beta\to -\infty}\log\frac{(\alpha-a)^2+b^2}{(\beta-a)^2+b^2}+ a(\tan^{-1}\alpha-\tan^{-1}\beta) \Big\}}}$$
$${\displaystyle{=\frac{1}{\pi}\Big\{ \lim_{\alpha\to\infty,\beta\to -\infty}\log\frac{(\alpha-a)^2+b^2}{(\beta-a)^2+b^2} \Big\} + \frac{1}{\pi}a\pi}}$$
これから、期待値が値を持たず、$${\alpha=-\beta}$$の時にのみ値$${a}$$を持つことがわかる。
コーシー分布は二次以上のモーメントは持たない。
$${\alpha=-\beta}$$の時の期待値の値を主値と呼ぶ。

標準コーシー分布$${a=0,b=1}$$に従う乱数の実装は、numpy.randomのstandard_cauchyで与えられる。

import numpy as np
from numpy import random as rand
import matplotlib.pyplot as plt
rand.seed(0)
randnum=1_000_000

cau=rand.standard_cauchy(size=randnum)
cau= cau[(cau>-5) & (cau<5)] 
plt.hist(cau, 'auto',label='standard caucy')

plt.ylabel('Probability')
plt.xlabel('Data')
plt.legend();

標準コーシー分布

任意の$${a,b\ge 0}$$のコーシー分布に従う乱数は、scipy.statのcaucyで実装することができる。

from scipy.stats import cauchy
mu=0.5
scale=1.0
cau=cauchy.rvs(loc=mu, scale=scale, size=randnum)
cau = cau[(cau>-5) & (cau<5)] 
plt.hist(cau, 'auto',label=f'$1/b$={scale}',alpha=0.6)
scale=2.0
cau=cauchy.rvs(loc=mu, scale=scale, size=randnum)
cau = cau[(cau>-5) & (cau<5)] 
plt.hist(cau, 'auto',label=f'$1/b$={scale}',alpha=0.6)
plt.title(f'Caucy distribution with a={mu}' )
plt.ylabel('Probability')
plt.xlabel('Data')
plt.legend();

ラプラス分布

母数$${1}$$の指数分布に従うそれぞれが独立な二つに確率変数の$${y,\ y'}$$の差$${\Delta y=y-y'}$$はは標準ラプラス分布に従い、その確率密度関数は
$${f(x)=\displaystyle{\frac{1}{2}exp(-|\Delta y|))}}$$と与えられる。
一般のラプラス分布の確率密度関数は$${f(x)=\displaystyle{\frac{1}{2b}exp(-\frac{|x-\mu|}{b})}}$$と与えられる。
モーメント母関数は、$${\displaystyle{|t|<\frac{1}{b}}}$$の時、
$${M_x(t)=\displaystyle{\frac{1}{2b}\int^a_{-\infty}exp(xt+\frac{x}{b}-\frac{\mu}{b})dx + \frac{1}{2b}\int^\infty_\mu exp(xt-\frac{x}{b}+\frac{\mu}{b})dx}}$$

$${=\displaystyle{\frac{e^{-\frac{\mu}{b}}}{2b}\Big[ \frac{b}{bt+1} exp(\frac{1+bt}{b}x)\Big]^\mu_{-\infty} +\frac{e^{\frac{\mu}{b}}}{2b}\Big[ \frac{b}{bt-1} exp(\frac{bt-1}{b}x)\Big]^{\infty}_\mu }}$$

$${=\displaystyle{\frac{1}{2} \Big(\frac{e^{\mu t}}{bt+1}-\frac{e^{\mu t}}{bt-1}\Big)=\frac{e^{\mu t}}{1-b^2t^2}}}$$
これから期待値は、
$${E[x]=\displaystyle{\frac{d}{dt}\frac{e^{\mu t}}{1-b^2t^2}\Big|_{t=0}=\Big(\frac{\mu e^{\mu t}}{1-b^2t^2}+\frac{2b^2te^{\mu t}}{(1-b^2t^2)^2}\Big)\Big|_{t=0} =\mu}}$$
分散は、
$${E[x^2]=\displaystyle{\frac{d^2}{dt^2}\frac{e^{\mu t}}{1-b^2t^2}\Big|_{t=0}}}$$

$${\displaystyle{=\Big( \frac{\mu^2 e^{\mu t}}{1-b^2t^2} +\frac{2b^2t\mu e^{\mu t}}{(1-b^2t^2)^2} + \frac{2b^2te^{\mu t} + 2\mu b^2t e^{\mu t}}{(1-b^2t^2)^2} +\frac{8b^2t^2e^{\mu t}}{(1-b^2t^2)^3}\Big)\Big|_{t=0} }}$$
$${=\mu^2+2b^2}$$から、$${V[x]=\mu^2+2b^2-\mu^2=2b^2}$$

これらのコーシー分布とラプラス分布は、異常値、外れ値を持つデータのモデル化に使われる。

b=1.0
mu=0.
lap=rand.laplace(loc=mu, scale=b, size=randnum)
plt.hist(lap, 100,label=f'b={b}',alpha=0.6)
b=2.0
lap=rand.laplace(loc=mu, scale=b, size=randnum)
plt.hist(lap, 100,label=f'b={b}',alpha=0.6)
b=4.0
lap=rand.laplace(loc=mu, scale=b, size=randnum)
plt.hist(lap, 100,label=f'b={b}',alpha=0.6)

plt.title(rf'Laplace distribution $\mu$={mu}')
plt.ylabel('Probability')
plt.xlabel('Data')
plt.legend();