PythonでRっぽいことができるpandas
Pythonにはpandasという素敵なライブラリがあります。pandasを簡単に言うと、PythonにRのデータフレームのような型を持たせるライブラリです。これによって、行列計算が劇的に楽になり、これまでRでやっていたような集計作業がPythonでも楽にできるようになります。
詳しくは公式ページをどうぞ
RでできるならRでやればいいじゃん、という意見はもっともですが、Pythonを使うことによるメリットも当然あります。
最も大きなメリットは、データクリーニングから集計・分析までPythonのみで実施できるという点です。データクリーニングに必要なテキスト処理などの地道な作業はPythonの方が得意です。当然、データベースからデータを抽出し、加工することもPythonなら普通にできます。もちろん、Rでもできることはできますが、Pythonほど便利ではありません。
地味なメリットとしては計算が早いという点です。ベンチマークをとったわけじゃないですが、体感的にはpandasで集計した方がかなり早いです。というのも、pandasでは裏でnumpyを利用しており、このnumpyはC言語で書かれています。また、公式サイトやオライリー本を読む限り、かなり最適化を行なっているようです。
計算処理が早いということは、古いPCでも大量のデータを扱えるということでもあります。SQLのJoinのようなデータの結合処理をした時は10万件以上のデータ同士でも、それほど時間がかからずに完了しました(ちなみに同じ作業をRで行ったら、PCがフリーズしました)。データ分析をする時は、速いが正義となります。
ただ、pandasにも弱点は当然あります。一番大きな問題はpandasには統計解析処理が「無い」という点です。。。集計まではできるのですが、統計処理を行うのは、scipyという別のライブラリでやれというスタンスです。scipyやその関連ライブラリも非常に充実しているのですが、Rと比べると弱かったり、ドキュメントが少なかったりして困ることが多いです。今回のお題のカイ2乗検定のやり方も中々見つからなくて困りました。
pandasでクロス集計
pandasの機能は膨大で全部紹介するのは無理ですので、今回はアンケート調査などでよくやるクロス集計表とその検定を行おうと思います。データはこちらです。男女別に好きな動物を聞きました的なイメージです。
データを読み込みます
import numpy as np
import pandas as pd
import scipy as sp
import scipy.stats
from pandas import DataFrame, Series
#データ読み込み
data = pd.read_csv("data.txt", sep="\t")
#データのサマリ
data.describe()
#=>
# gender animal
#count 38 38
#unique 2 3
#top f dog
#freq 25 15
#単純集計
for col in data:
print pd.value_counts(data[col])
print "\n"
#=>
#f 25
#m 13
#dtype: int64
#
#dog 15
#cat 14
#usagi 9
#dtype: int64
こんな感じで簡単にデータを読み込んで、サマリを見ることができます。
pandasでクロス集計
クロス集計も簡単です。ついでに行%表も作っちゃいましょう。
crossed = pd.crosstab(data.gender, data.animal)
print crossed
#=>
#gender cat dog usagi
#f 3 13 9
#m 6 2 5
#行%表の作成
arr = []
for row in crossed.T:
arr.append(crossed.T[row] / float(crossed.T[row].sum()))
crossed_per = DataFrame(arr)
print crossed_per
#=>
#animal cat dog usagi
#f 0.120000 0.520000 0.360000
#m 0.461538 0.153846 0.384615
ここまではpandasのみでサクッと作れちゃいます。
このクロス集計表に対して、カイ2乗検定を行います。
scipyでカイ2乗検定
クロス集計表のカイ2乗検定はscipy.stats.chi2_contingencyを使います。
scipy.stats.chi2_contingencyの戻り値は、tuple型で上から
・カイ2乗値
・P値
・自由度
・期待値
となっています。
ではカイ2乗検定を行います。
#p-value
print res[1]
#=> 0.028280095
P値が0.05以下なので、有意差ありという結果です。
地道に残差分析
ここからが本番です。
カイ2乗検定で有意差が認められた場合、残差分析を行なってクロス表のどのセルが効いて違いを作ったのかを探ります。
Rであればchisq.test()の戻り値に調整済み標準化残差が含まれているので簡単なのですが、scipyには探した限りありません。
自分で作るしかありません。。。
res = {}
#カイ2乗検定の実施→カイ2乗値、p値、自由度、期待値が戻り値
df_chi = sp.stats.chi2_contingency(df)
res = {
"data" : df,
#p値
"p_val" : df_chi[1],
#期待値
"df_exp": DataFrame(df_chi[3])
}
#期待値のカラム名とインデックス名を基データに合わせる
res["df_exp"].columns = df.columns
res["df_exp"].index = df.index
#残差
res["df_res"] = df - res["df_exp"]
res["df_res"].columns = df.columns
res["df_res"].index = df.index
#行%の計算
arr = []
for row in df.T:
arr.append(df.T[row] / float(df.T[row].sum()))
res["df_per"] = DataFrame(arr)
res["df_per"].columns = df.columns
res["df_per"].index = df.index
#残差分析用前処理
row_sum = df.T.sum()
col_sum = df.sum()
full_sum = float(row_sum.sum())
#残差分散を算出
arr_all = []
for r in row_sum:
arr = []
for c in col_sum:
arr.append((1-(r/full_sum))*(1-(c/full_sum)))
arr_all.append(arr)
res["df_res_var"] = DataFrame(arr_all)
res["df_res_var"].columns = df.columns
res["df_res_var"].index = df.index
col_size = df.columns.size
row_size = df.index.size
#調整済み標準化残差を算出
arr_all = []
for r in np.arange(row_size):
arr = []
for c in np.arange(col_size):
arr.append(res["df_res"].iloc[r].iloc[c] / np.sqrt(res["df_exp"].iloc[r].iloc[c] * res["df_res_var"].iloc[r].iloc[c]))
arr_all.append(arr)
res["df_res_final"] = DataFrame(arr_all)
res["df_res_final"].columns = df.columns
res["df_res_final"].index = df.index
return res
#関数の実行
res = chi_sq_test(crossed)
#p値の出力
print(res["p_val"])
#=> 0.028280095
#調整済み標準化残差の出力
print(res["df_res_final"])
#=>
#gender cat dog usagi
#f -2.349377704162738 2.190723321954562 -0.14923492309463762
#m 2.349377704162738 -2.1907233219545614 0.149234923094637
catとdogが1.96を超えています。この2つが影響を与えて有意差が出たと言えそうです。
あまりPythonに慣れていないので、残差分散とか調整済み標準化残差とかの計算がもっとエレガントにできそうな気はしますが、一応これで残差分析までできました。
データ集計までは非常に楽できるのですが、そこからscipy使うのに大きなハードルがあるという感じですね。単純にインタラクティブな分析をするのであれば、Rの方が使いやすいですが、アプリケーションとして実行させるならPythonのみで作った方が楽そうです。
5 replies on “Pythonデータ解析ライブラリpandasと遊ぶ:クロス集計~検定・残差分析まで”
[…] こちらのサイト Pythonデータ解析ライブラリpandasと遊ぶ:クロス集計~検定・残差分析まで のCSVファイルを読み込んでみる。 […]
h Black Dial Two Tone
Fairy Fashion D
phone 6(4.7) with Screen
[…] http://web-analytics-or-die.org/2013/07/pandas/ […]