KaggleでEDAをやろうと思った時に見る記事

#記事を書いた目的
Kaggleコンペに参加する時、まず最初にやるべしと言われるのがEDA(Explanatory Data Analysis＝探索的データ分析)です。ただ、何をやれば良いのかあまりまとまったものがないと思ったので、備忘のためにもまとめてみました。PythonとJupyter Notebookを使って数値予測のテーブルデータコンペへ参加する前提です。

こちらのKernelを参考にしました。

間違っている点や追加した方が良い点など、アドバイスをコメントで頂けると大変ありがたいです。

#EDAでやること
・データのカラム名と意味の確認
・データの大きさとデータ型の確認
・欠損値チェック
・カラムごとの基礎統計量確認
・ヒストグラムによる分布の可視化
・散布図や箱ひげ図を用いた目的変数と、説明変数の関係の可視化
・時系列データであれば、時系列トレンドの可視化
・外れ値の確認
あたりでしょうか。

#必要なライブラリのインポート

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
from scipy import stats
%matplotlib inline #matplotlibの図をJupyter Notebook上で表示するためのマジックコマンド

#データの読み込み
pd.read_csv("ファイルパス")を使います、KaggleのKernelだと大体以下で上手くいくはず。

train = pd.read_csv("../input/train.csv")
test = pd.read_csv("../input/test.csv")

#データのカラム名と意味の確認
KaggleだとCompetitionのDataタブで詳細情報が出てくるので、これを熟読します。
これと併せて、OverviewにあるEvaluationも見ておくことをお勧めします。どのようなMetricで評価されるか理解しておくことはコンペへ参加する上で非常に大切です。

#データの大きさとデータ型の確認
データフレーム.shapeとデータフレーム.info()で大体いけます。
info()では欠損値でないものの数も確認できるため、結果的に欠損値の多いカラムも分かります。

print(train.shape)
print("*"*50)
print(test.shape)
print("*"*50)
print(train.info())
print("*"*50)
print(test.info())

#データの確認
データフレーム.head()でデータそのものも見ておきましょう。
Jupyter Notebookでは、結果を表示するためにセルを分けて実行しましょう。

train.head()

test.head()

また、Numerical featuresとCategorical featuresのカラム名をそれぞれ持っておくと、何かと便利です。

numerical_feats = train.dtypes[train.dtypes != "object"].index
categorical_feats = train.dtypes[train.dtypes == "object"].index

例えばtrain[numerical_feats]とすればNumerical featuresだけを表示することが出来ます。

#欠損値チェック
参考先のコードがそのまま使えるかと。

total = train.isnull().sum().sort_values(ascending=False)
percent = (train.isnull().sum()/train.isnull().count()).sort_values(ascending=False)
missing_data = pd.concat([total, percent], axis=1, keys=['Total', 'Percent'])
missing_data.head(20)

データフレーム.isnull(): 欠損値であればTrue、それ以外はFalseを入れたデータフレームを返す。
データフレーム.sum() : 数値の合計を求める。Boolean型だとTrueは1、Falseは0となる。
データフレーム.count() : データ数を求める。
データフレーム.sort_values() : データを並び替える。ascending=Falseだと降順となる。
pd.concat([データフレーム1, データフレーム2], axis=1, keys=['カラム名1', 'カラム名2']) : axis=1とすることで2つのデータフレームを横方向で結合する。keys以下で新たなカラム名を振っている。

これはつまり、percentの分母であるtrain.isnull().count()はデータ行数のため全カラム同一であることから、totalの降順＝percentの降順となり、2つのデータフレームを横方向で結合することが出来るということです。

なお、EDAの次のデータ前処理の段階で、欠損値を何らかの値で置き換えたり、欠損値の多いカラムを削除したりすることが必要となります。

train[col].fillna(NaNを置き換える値, inplace=True)
#inplace=Trueで元のデータフレーム自体を書き換える

のように書きますが、データフレーム中のNaNが意味のある情報(参考先だとno poolingを意味する等)なのか、本当の欠損値なのかは、注意する必要がある模様です。

#Numerical features編
##基礎統計量確認
データフレーム.describe()で一発です。
こちらもJupyter Notebookでは、結果を表示するためにセルを分けて実行しましょう。

train.describe()

test.describe()

##ヒストグラムによる分布の可視化
カラム毎に

sns.distplot(train['カラム名'])

としてヒストグラムを見ていくことになります。
参考先にあるように

for col in numerical_feats:
    print('{:15}'.format(col), 
          'Mean: {:05.2f}'.format(train[col].mean()) , 
          '   ' ,
          'Std: {:05.2f}'.format(train[col].std()) , 
          '   ' ,
          'Skewness: {:05.2f}'.format(train[col].skew()) , 
          '   ' ,
          'Kurtosis: {:06.2f}'.format(train[col].kurt())  
         )

としてカラム毎の平均、標準偏差、歪度、尖度を求めるのも参考になるかと思います。
正規分布っぽくないやつがあれば対数変換することも選択肢です。

##目的変数と説明変数の関係性と外れ値の確認
参考先の以下コードがとても良さ気です。

nr_rows = 12 #図を表示する際の行数
nr_cols = 3 #図を表示する際の列数
#nr_rows * nr_colsがカラム数を超えるように設定。

fig, axs = plt.subplots(nr_rows, nr_cols, figsize=(nr_cols*3.5,nr_rows*3))

li_num_feats = list(numerical_feats)
li_not_plot = [説明変数として使わないカラム名(IDや目的変数など)]
li_plot_num_feats = [c for c in list(numerical_feats) if c not in li_not_plot]
target = "目的変数のカラム名"

for r in range(0,nr_rows):
    for c in range(0,nr_cols):  
        i = r*nr_cols+c
        if i < len(li_plot_num_feats):
            sns.regplot(train[li_plot_num_feats[i]], train[target], ax = axs[r][c])
            stp = stats.pearsonr(train[li_plot_num_feats[i]], train[target])
            str_title = "r = " + "{0:.2f}".format(stp[0]) + "      " "p = " + "{0:.2f}".format(stp[1])
            axs[r][c].set_title(str_title,fontsize=11)
            
plt.tight_layout()    
plt.show()  

R値が大きいものは説明変数として使える可能性が高いと言えるでしょう。
また、各散布図を眺めることで外れ値を見つけることが出来ます。
外れ値となるデータ行の削除は

train = train.drop(train[外れ値の条件].index)

とします。

#Categorical features編
##基礎統計量確認
基礎統計量としてUnique valueを見ます。参考先のコードがそのまま使えるかと。

for catg in list(categorical_feats) :
    print(train[catg].value_counts())
    print('#'*50)

##目的変数と説明変数の関係性
参考先の以下コードがとても良さ気です、箱ひげ図をまとめて描きます。

li_cat_feats = list(categorical_feats)
nr_rows = 15 #図を表示する際の行数
nr_cols = 3 #図を表示する際の列数
#nr_rows * nr_colsがカラム数を超えるように設定。

fig, axs = plt.subplots(nr_rows, nr_cols, figsize=(nr_cols*4,nr_rows*3))

for r in range(0,nr_rows):
    for c in range(0,nr_cols):  
        i = r*nr_cols+c
        if i < len(li_cat_feats):
            sns.boxplot(x=li_cat_feats[i], y=target, data=train, ax = axs[r][c])
    
plt.tight_layout()    
plt.show()  

あるカラムにおいてUnique value毎に中央値や箱の範囲が大きくずれているほど、各Unique valueとTargetの関係性が強いため、予測に使えるカラムだと言える。逆に、どのUnique valueでも箱が同じような範囲になっていると、予測に使えない。

#その他
・この後は、使えるカラムを抜き出しつつ、データ前処理を行い、モデル構築へつなげていく。
・使えるカラムの抽出では、Numerical featuresでは相関係数で閾値を設ける、Categorical featuresではまずは箱ひげ図の目視、で良い？
・Categorical featuresを予測で用いるには数値へのencodingが必要。単純なLabel encodingでもTree based model(XgboostやLightGBMなど)では使えるが、Unique valueが出現する頻度にデータを置き換えるFrequency encodingも役立つケースがあるし、Tree based model以外ではOne-Hot Encodingが必要なケースが多い。
・pandas-profilingが便利らしい。
・Categorical featuresについてはこのKernelも参考になりそう、Violin plotの使い方なども別途まとめてみたい。

#(参考)データフレームのメモリ使用量圧縮(reduce_mem_usage)
Kaggleの有志が作成した関数が便利なので、こちらにあったコードを書いておきます。

def reduce_mem_usage(props):
    start_mem_usg = props.memory_usage().sum() / 1024**2 
    print("Memory usage of properties dataframe is :",start_mem_usg," MB")
    NAlist = [] # Keeps track of columns that have missing values filled in. 
    for col in props.columns:
        if props[col].dtype != object:  # Exclude strings
            
            # Print current column type
            print("******************************")
            print("Column: ",col)
            print("dtype before: ",props[col].dtype)
            
            # make variables for Int, max and min
            IsInt = False
            mx = props[col].max()
            mn = props[col].min()
            
            # Integer does not support NA, therefore, NA needs to be filled
            if not np.isfinite(props[col]).all(): 
                NAlist.append(col)
                props[col].fillna(mn-1,inplace=True)  
                   
            # test if column can be converted to an integer
            asint = props[col].fillna(0).astype(np.int64)
            result = (props[col] - asint)
            result = result.sum()
            if result > -0.01 and result < 0.01:
                IsInt = True

            
            # Make Integer/unsigned Integer datatypes
            if IsInt:
                if mn >= 0:
                    if mx < 255:
                        props[col] = props[col].astype(np.uint8)
                    elif mx < 65535:
                        props[col] = props[col].astype(np.uint16)
                    elif mx < 4294967295:
                        props[col] = props[col].astype(np.uint32)
                    else:
                        props[col] = props[col].astype(np.uint64)
                else:
                    if mn > np.iinfo(np.int8).min and mx < np.iinfo(np.int8).max:
                        props[col] = props[col].astype(np.int8)
                    elif mn > np.iinfo(np.int16).min and mx < np.iinfo(np.int16).max:
                        props[col] = props[col].astype(np.int16)
                    elif mn > np.iinfo(np.int32).min and mx < np.iinfo(np.int32).max:
                        props[col] = props[col].astype(np.int32)
                    elif mn > np.iinfo(np.int64).min and mx < np.iinfo(np.int64).max:
                        props[col] = props[col].astype(np.int64)    
            
            # Make float datatypes 32 bit
            else:
                props[col] = props[col].astype(np.float32)
            
            # Print new column type
            print("dtype after: ",props[col].dtype)
            print("******************************")
    
    # Print final result
    print("___MEMORY USAGE AFTER COMPLETION:___")
    mem_usg = props.memory_usage().sum() / 1024**2 
    print("Memory usage is: ",mem_usg," MB")
    print("This is ",100*mem_usg/start_mem_usg,"% of the initial size")
    return props, NAlist

reduce_mem_usage(データフレーム)として実行すれば、Numerical featuresのデータ型を適切なものへ変換することで、メモリ使用量を圧縮してくれます。データサイズの大きいコンペで活躍します。