Search on the blog

2015年12月2日水曜日

カイ二乗値を計算する

 カイ二乗値を計算して、カイ二乗検定とはそういうことだったのかと今さら納得した。

まず、偏りのないサイコロを振ってカイ二乗値を計算してみる。
import bisect

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import chi2

np.random.seed(1226)

# create a dice with given frequency
def create_dice(freq):
    cdf = np.cumsum(freq)
    def dice():
        p = np.random.uniform()
        return bisect.bisect_left(cdf, p)
    return dice

# calculate chi-squared value
def calc_chi_sq(freq, sample_size):
    dice = create_dice(freq)
    sample = np.bincount([dice() for _ in xrange(sample_size)])
    populatoin = [sample_size * f for f in freq]
    return np.sum([(x-y)**2/y for (x, y) in zip(sample, populatoin)])

freq = [1./6, 1./6, 1./6, 1./6, 1./6, 1./6]

# plot dice histogram
dice = create_dice(freq)
plt.hist([dice() for _ in range(3000)], bins=range(7), normed=True)
plt.show()

# plot chi-squared
chi_sqs = [calc_chi_sq(freq, 3000) for _ in xrange(3000)]
plt.hist(chi_sqs, bins=50, normed=True)
x = np.linspace(0,35,50)
y = chi2.pdf(x, df=5)
plt.plot(x, y, color='r')
plt.show()

実行したら以下の2つのグラフが出力される。
上がサイコロの目の分布。下がカイ二乗値の分布。


カイ二乗値はカイ二乗分布に従っていることが分かる。

続いて、激しく偏った分布を持つサイコロを振って、同じことをしてみる。
上のコードの


freq = [1./6, 1./6, 1./6, 1./6, 1./6, 1./6]
の部分を
freq = [3./12, 1./12, 4./12, 2./12, 1./12, 1./12]
に変えて実行してみると、


なるほど、母集団の分布に関係なくカイ二乗値はカイ二乗分布に従うと。
ということで、カイ二乗検定が何をやってるのか今さらながら納得できた。。。

0 件のコメント:

コメントを投稿