スポンサーリンク
スポンサーリンク

カイ二乗分布が自然に理解できる記事

統計学 統計学

機械学習を先に学んでしまうと、カイ二乗分布が良く出てくる理由が分からないかもしれません。そこで、機械学習(ベイズ統計)では良く使うガンマ分布の一つとして定式化しておいて、カイ二乗分布の大事さ、特に標準正規分布との関係を見ていきます。
参考文献は一応2冊あります。

統計学入門 (基礎統計学Ⅰ) | 東京大学教養学部統計学教室 |本 | 通販 | Amazon
Amazonで東京大学教養学部統計学教室の統計学入門 (基礎統計学Ⅰ)。アマゾンならポイント還元本が多数。東京大学教養学部統計学教室作品ほか、お急ぎ便対象商品は当日お届けも可能。また統計学入門 (基礎統計学Ⅰ)もアマゾン配送商品なら通常配送無料。
https://amzn.to/2WjGIbY
スポンサーリンク

モーメント母関数と特性関数

初めに、話をスムーズにするためにモーメント母関数と特性関数を定義しておきます。
Xが確率変数で、どの確率分布を考えているか明らかな時、Xの期待値をE[X]で表します。モーメント母関数とは、次の量を指します。
MX(t)=E[exp(Xt)]
モーメント母関数の名前の由来は、モーメント母関数から全ての次数のモーメントが得られる事から来ています。1
例えば、正規分布N(μ,σ2)のモーメント母関数は、
MX(t)=exp(μt+σ2t2/2)
です。一般に、モーメント母関数は存在するとは限りませんが、以下で定義する、特性関数はいつも存在します。
ϕX(t)=E[exp(iXt)]
フーリエ変換、逆変換の関係から、特性関数と確率分布は1対1に対応します。また、存在すれば、特性関数とモーメント母関数の間には以下の関係があるので、モーメント母関数からも確率分布は確定します。
ϕX(t)=MX(it)
このことから、確率分布の性質を議論する上で、楽が出来る事があります。例えば、正規分布では
N(μ1,σ12)+N(μ2,σ22)=N(μ1+μ2,σ12+σ22)
が成り立ちますが、モーメント母関数を使うと楽になります。
独立な確率変数をXN(μ1,σ12),YN(μ2,σ22)として、
Z=X+YN(μ1+μ2,σ12+σ22)が成り立つ事を確かめます。
MZ(t)=E[exp(Zt)]=E[exp(Xt)]E[exp(Yt)]=exp(μ1t+σ12t2/2)exp(μ2t+σ22t2/2)=exp((μ1+μ2)t+(σ12+σ22)t2/2)
これは、正規分布 N(μ1+μ2,σ12+σ22) のモーメント母関数です。

ガンマ分布とカイ二乗分布

ガンマ分布とカイ二乗分布について簡単な性質を計算してみます。ガンマ分布は、以下のような分布でした。色々な計算をまとめた記事があるので、計算はそちらを見てください。
Gam(x|α,β)=1Γ(α)1β(xβ)α1ex/β

ガンマ分布たち


ガンマ分布のモーメントと母関数や、期待値、分散を書いておきます。
Mx(t)=(1βt)αE[X]=αβVar[X]=αβ2
モーメント母関数の形から分かるように、独立なガンマ分布にも、ある種の再生性があります。
Gam(α1,β)+Gam(α2,β)=Gam(α1+α2,β)
一応計算しておきましょう。XGam(α1,β),YGam(α2,β)で、XとYは独立とします。
E[exp((X+Y)t)]=E[exp(Xt)]E[exp(Yt)]=(1βt)α1(1βt)α2=(1βt)(α1+α2)
最後の式は、
Gam(α1+α2,β)のモーメント母関数なので、欲しい結論が得られました。

ガンマ分布の特別な場合をカイ二乗分布と呼びます。

Gam(x|n/2,2)を自由度n のカイ二乗分布と呼び、χn2で表す。

カイ二乗分布
カイ二乗分布

カイ二乗分布の確率分布や、モーメント母関数などを表すと、以下のようになります。
χn2 =1Γ(n/2)12(x2)n/21ex/2Mx(t)=(12t)n/2E[X]=nVar[X]=2n
ガンマ分布なので、カイ二乗分布の間でも再生性があります。
χn2+χm2=χn+m2
わざわざ名前が付いている分布なので、大事な分布なのですが、そのことを見ていきます。

標準正規分布とカイ二乗分布

カイ二乗分布が登場するシチュエーションを見てみます。

ZN(0,1)のとき、Z2χ12に従う。

標準正規分布があれば、カイ二乗分布があるのです。一応計算しておきます。
モーメント母関数を計算して、χ12が出てくることを確かめます。
MZ2(t)=E[ez2t]=12πdzez2tez2/2dz=12πdzexp(12t2z2)dz=12π112tdzexp(x22)dx=(12t)1/2
途中で、x=(12t)zと変数変換しました。
最後の式はχ12のモーメント母関数なので、Z2が自由度1のカイ二乗分布に従う事が分かりました。
上の主張と、自由度の足し算が出来るということから、次の事が分かります。

Z1,,ZnN(0,1)がそれぞれ独立なとき、Zi2χn2に従う。

実際に標準正規分布からのサンプルを二乗して足したものをヒストグラムに描いてみると、カイ二乗分布とほぼ一致することが確かめられます。

np.random.seed(0)
z=np.zeros([10,100])
for i in range(10):
    z[i]=np.random.randn(100)
chi_ = np.sum(z**2, axis=0)

chi = np.random.chisquare(10,100)

sns.distplot(chi_,label="10個の標準正規分布から")
sns.distplot(chi,label="自由度10カイ二乗分布から")
plt.title("分布の比較")
plt.legend()
plt.savefig("chi_from_norms.png")
正規分布からカイ二乗分布
正規分布からカイ二乗分布

さらに大事な事として、正規分布からサンプリングして、不偏分散を考えるとカイ二乗分布が出てくるという話があります。

不偏分散とカイ二乗分布

X1,,XnN(μ,σ2)とする。また、それぞれの変数は独立とする。Zi=(XiX¯)/σとする。ただし、X¯Xi 達の標本平均で、X¯=Xi/nを表す。
このとき、
i=1n(ZiZ¯)2χn12
が成り立つ。2
さらに、 Z¯i=1n(ZiZ¯)2は独立。

上の主張のキモは、自由度がnでなくて、n-1 になるという事です。
i=1n(ZiZ¯)2=Zi2nZ¯2であり、Zi,nZ¯は標準正規分布に従いますが、nZ¯2Zi2達の自由度を奪っています。これは証明が面倒なので、大体の流れだけを書いておきます。
正確な証明は、参考文献を読んでください。

[証明の概略]
Ziなどが標準正規分布になる事は、モーメント母関数を計算すれば分かる。
Z1,,Zn,Z¯の組で考えると、どれが独立かはっきりしないので、直交変換Hで、(Z1,,Zn)(Y1=Z¯,Y2,,Yn)となるようなものを作る。3
これが作れると、
i=1n(ZiZ¯)2=Zi2nZ¯2=ZTZY12=(HY)THYY12=i=2nYi2
となる。Zi同士は独立で、標準正規分布に従い、Hが直交変換なので、Yiも標準正規分布に従う。
最後に、Y1,,Ynは全て独立で、 上の計算から、
i=2nYi2=i=1n(ZiZ¯)2Y1=nZ^
なので、Z¯i=1n(ZiZ¯)2は独立。

実際に、統計量
(n1)Var(X)/σ2
を作ってヒストグラムにしてみると、カイ二乗分布とほぼ一致することが確認できます。

#データの準備 10種類の正規分布から100組のデータを集める
sample=10
size=100
sig=2
mu=1
x=np.zeros((size,sample))
for i in range(sample):
    np.random.seed(i)  
    x[:,i]=np.random.normal(loc=mu,scale=sig,size=size)
#統計量の計算とカイ二乗分布の準備
chi_=[]
for j in range(size):
    v=np.sum( (x[j,:] - np.mean(x[j,:]))**2)/(sample-1)
    V = (sample-1)*v/sig**2
    chi_.append(V)

chi = np.random.chisquare(9,100)
#グラフを描く
sns.distplot(chi_,label="10個の正規分布から")
sns.distplot(chi,label="自由度9カイ二乗分布から")
plt.title("分布の比較")
plt.legend()
plt.savefig("chi_from_norms2.png")
不偏分散からカイ二乗分布が作られる
不偏分散からカイ二乗分布が作られる

このような背景があるので、色々な場所でカイ二乗分布が出てくるわけです。
しかし、いきなり正規分布が得られている事は少ないので、近似的に正規分布が得られて、その結果カイ二乗分布が出てくる、という状況の方が多いかもしれません。
そのような状況の記事は別に書きたいと思います。

多変量正規分布とカイ二乗分布

標準正規分布に従う確率変数X1,,XnN(0,1)があるとき、自由度nのカイ二乗分布が出て来ました。
i=1nXinχn2
これらを、確率変数からなるベクトルだと思うと、以下のように書き直すことが出来ます。
X=(X1,,Xn)XTX=i=1nXinχn2
多変量正規分布の話では、確率変数ベクトルの成分の共分散が0なら、成分は独立です。
つまり、XN(0,In)なら、
XiN(0,1)XTX=i=1nXi2χn2
という事です。
一変数の場合のように、標準多変量正規分布からカイ二乗分布が出て来ます。

X=(X1,,Xn)N(μ,Σ) とします。Xを標準化すると、自由度nのカイ二乗分布が得られます。4
(Xμ)TΣ1(Xμ)χn2

これを証明するには、以下の事実を使います。5
XN(μ,Σ)AX+bN(Aμ+b,AΣAT)
上の事実を認めると、
XμN(0,Σ)
です。なので、XN(0,Σ) に対して、
AΣAT=IXTATAX=XTΣ1X
となるようなAが存在する事を示せば、AXN(0,AΣAT)=N(0,I)となり、主張が示せます。
正規分布に使われる共分散行列が、正定値かつ対称としましょう。6
そうすると、Σは直交行列で対角化できます。
HΣHT=ΛHTH=I
ただし、Λ=diag[λ1,,λnは固有値が対角に並んだ行列です。Λ=diag[1/λ1,,1/λn]とおいて、A=HΛとしましょう。
AΣAT=Λ2Λ=I
であり、Σ=HTΛHなので、ATA=Σ1が分かります。この事から
XTATAX=XTΣ1X
となります。これで、主張が示せました。
python で主張を確かめてみます。

#データの生成
np.random.seed(2)
sig=np.random.randint(1,10,5)
Sigma = np.dot(sig.reshape(-1,1), sig.reshape(1,-1)) + 1*np.eye(5)
mu =np.ones(5)
y=np.random.multivariate_normal(mean=mu ,cov=Sigma, size=100)
#
Y=np.array([])
for i in range(len(y)):
    chi_ = (y-mu)[i].reshape(-1,5)@ np.linalg.inv(Sigma)@(y-mu)[i].reshape(5,-1)
    Y=np.append(Y,chi_ )

chi = np.random.chisquare(5,100)
#グラフを描く
sns.distplot(Y,label="5次元正規分布から")
sns.distplot(chi,label="自由度5カイ二乗分布から")
plt.title("分布の比較")
plt.legend()
plt.savefig("chi_from_norm.png")
5次元正規分布→カイ二乗分布
5次元正規分布から自由度5のカイ二乗分布が出てくる

上の主張
(Xμ)TΣ1(Xμ)χn2
は、変量正規分布を標準化するとカイ二乗分布が出てくるという意味ですが、一方で正規分布の肩に載っている量がカイ二乗分布になるという事でもあります。
(Xμ)TΣ1(Xμ)は、Xが正規分布N(μ,Σ)からの距離と捉える事も出来ます。7
この量を使って、異常検知のシステムを作る事が出来ます。
このことは別の記事で解説します。

まとめ

  • 計算するためにモーメント母関数などを定義した
  • ガンマ分布の一つとしてカイ二乗分布を定義した
  • 標準正規分布からカイ二乗分布を作った
  • 分散からカイ二乗分布を作った
  • 多変量正規分布からカイ二乗分布を作った。
  1. E[X],E[X2],,を、1次,2次,のモーメントと呼びますが、モーメント母関数をX=0 でn階微分すれば、n次モーメントが得られます。
  2. (n1)Var(X)/σ2χn12という事です。
  3. Z=HY かつ、HTH=IとなるHを作ります。Helmert 行列というらしいです。
  4. nは、ベクトルの成分の数です。
  5. 多変量正規分布の記事で証明してます。
  6. 正定値は、全ての固有値が0より大きいという事です。
  7. マハラノビス距離と呼ばれます。
タイトルとURLをコピーしました