[Python]とにかくわかりやすく!ハンバーガー統計をやってみる②

14786 ワード
scipy Python 初心者 Python3 統計学 Python テキストリンク

前回の記事

前回の記事→こちら

環境

Mac OSX　mojave
Python３.7

とりあえず超超有名な入門書

ハンバーガー統計をPythonで書いてみたいと思います。
元はこちらの内容

分量が多くなるので
①：平均と分散、信頼区間
②：カイ二乗検定 ←本記事はここ
③：t検定
④：分散分析
のようなシリーズで行きたいと思います。

カイ二乗検定

ワクワクバーガーのチキンが、ポテトよりも売れていないくて、それがモグモグバーガーよりも売れていないかもという話から売上の比較を行うことになりました。

csvでデータを用意しました。

さて検定なので、原則は「差はない」という帰無仮説から始まります。

では計算をしていきます。
実際に売れた個数を観測度数と呼び、全体と同じ割合で売れるときの個数を期待度数と言います。期待度数を計算してみます。

まずはpivot_table()を使います。crosstabでもいいです。

uriage_dframe = pd.read_csv("uriage_data.csv")

cross = pd.pivot_table(
    data = uriage_dframe,
    values = "cnt",
    aggfunc = "sum",
    index = "store",
    columns = "item",
    margins = True)

ピボットテーブルができました。

それでは全体のポテトの割合を算出します。四則演算です。

#全体に置けるチキンの割合。ポテトの割合
potates_share = cross["potates"][2]/cross["All"][2]
chickens_share = cross["chickens"][2]/cross["All"][2]

この割合を元に、期待度数を算出します。

#全体と同じ割合で売れるとした時の個数
ex_freq_wakupotates = cross["All"][1]*potates_share
ex_freq_wakuchickens = cross["All"][1]*chickens_share
ex_freq_mogupotates = cross["All"][0]*potates_share
ex_freq_moguchickens = cross["All"][0]*chickens_share
print(ex_freq_wakupotates)
print(ex_freq_wakuchickens)
print(ex_freq_mogupotates)
print(ex_freq_moguchickens)
420.0
180.0
280.0
120.0

例えば、この観測度数４３５に対して期待度数420という、このズレが誤差なのかどうかをみてみようと思います。

カイ二乗値（ズレの大きさ）を算出します。
・$(観測度数-期待度数)$の総和　←０になる
・$((観測度数-期待度数)^2)$の総和　←考え方はOKだがデータ量が大きければ大きいほどズレが大きくなる
・$(((観測度数-期待度数)^2)\div{期待度数})$の総和 ←これがカイ二乗値

計算してみます。

a = (cross["potates"][1]-ex_freq_wakupotates)**2/ex_freq_wakupotates
b = (cross["potates"][0]-ex_freq_mogupotates)**2/ex_freq_mogupotates
c = (cross["chickens"][1]-ex_freq_wakuchickens)**2/ex_freq_wakuchickens
d = (cross["chickens"][0]-ex_freq_moguchickens)**2/ex_freq_moguchickens

chi_2 = a+b+c+d
print(chi_2)
4.464285714285714

カイ二乗分布は、大きくなればなるほど、確率密度は小さくなります(=起きにくくなります)
また自由度によって分布の仕方も変わります。分布表は別途確認してください。
$自由度 = (行の数−1)\times(列の数−1)$で求められます。

ここまで踏まえて計算してみます。

#自由度を求めます
df = (len(ex_freq_dframe.index)-1)*(len(ex_freq_dframe.columns)-1)

if chi_2>3.84:
    print("自由度が"+str(df)+"で有意水準を5%とした時に、カイ二乗値が"+str(chi_2)+"なので、帰無仮説を棄却し対立仮説を採択します")
else:
    print("自由度が"+str(df)+"で有意水準を5%とした時に、カイ二乗値が"+str(chi_2)+"なので、帰無仮説を採択します")

結果：自由度が1で有意水準を5%とした時に、カイ二乗値が4.464285714285714なので、帰無仮説を棄却し対立仮説を採択します

ということで有意水準（めったに起きない水準＝誤差とはいえない水準）よりも大きく　"「差がない」とは言えない"　という結果が出ました。

ちなみにカイ二乗値が出ていれば以下のようにもだせます。p値で処理しています。

p = 1 - sp.stats.chi2.cdf(x = chi_2 ,df = df)
if p < 0.05:
    print("結果：自由度が"+str(df)+"で有意水準を5%とした時に、カイ二乗値が"+str(chi_2)+"なので、帰無仮説を棄却し対立仮説を採択します")
else:
    print("結果：自由度が"+str(df)+"で有意水準を5%とした時に、カイ二乗値が"+str(chi_2)+"なので、帰無仮説を採択します")

結構コード長くなりましたが、これ一発で出ます。※事前にpivotのmargins=Falseにしないと自由度がおかしくなります。

cross = pd.pivot_table(
    data = uriage_dframe,
    values = "cnt",
    aggfunc = "sum",
    index = "store",
    columns = "item",
    margins = False)

result_chi2 = sp.stats.chi2_contingency(cross, correction = False)
print(result_chi2)
(4.464285714285714, 0.03461055751570723, 1, array([[120., 280.],
       [180., 420.]]))

カイ二乗値、p値、自由度、期待度数のarrayの順番で結果が出ています。

終わり

なんかすごく短い投稿になりました。続きはこちら