はじめに

Kaggle に少し手をつけてみようと思い EDA kernel (データを並べてみて解析するカーネル？) を眺めてみたのですが何をしているのかよくわからない部分が多かったので、ブログでまとめてみようと思います。
とはいえそれなりに量が多いので全てのコードは扱いませんし、一度に書き切れる気もしないのでちょくちょく更新しようと考えています。

参考にしたカーネルはこちらです

%%time
root = '../input/ashrae-energy-prediction/'
train_df = pd.read_csv(root + 'train.csv')
train_df["timestamp"] = pd.to_datetime(train_df["timestamp"], format='%Y-%m-%d %H:%M:%S')

weather_train_df = pd.read_csv(root + 'weather_train.csv')
test_df = pd.read_csv(root + 'test.csv')
weather_test_df = pd.read_csv(root + 'weather_test.csv')
building_meta_df = pd.read_csv(root + 'building_metadata.csv')
sample_submission = pd.read_csv(root + 'sample_submission.csv')

root で入力csv のルートを指定して、 pd.read_csv でcsvを読み込んでいるようですね。 pd.to_datetime でパースしているようですが、これは

note.nkmk.me

このサイトによると、

pandas.to_datetime()関数を使うと、日時（日付・時間）を表した文字列の列pandas.Seriesをdatetime64[ns]型に変換できる。と言うことだそうですね。 Panda.Series はある一列のことのよう。

あとは同じような感じで色々なcsvを読み込んでいると。

## Function to reduce the DF size
def reduce_mem_usage(df, verbose=True):
    numerics = ['int16', 'int32', 'int64', 'float16', 'float32', 'float64']
    start_mem = df.memory_usage().sum() / 1024**2    
    for col in df.columns:
        col_type = df[col].dtypes
        if col_type in numerics:
            c_min = df[col].min()
            c_max = df[col].max()
            if str(col_type)[:3] == 'int':
                if c_min > np.iinfo(np.int8).min and c_max < np.iinfo(np.int8).max:
                    df[col] = df[col].astype(np.int8)
                elif c_min > np.iinfo(np.int16).min and c_max < npfor.iinfo(np.int16).max:
                    df[col] = df[col].astype(np.int16)
                elif c_min > np.iinfo(np.int32).min and c_max < np.iinfo(np.int32).max:
                    df[col] = df[col].astype(np.int32)
                elif c_min > np.iinfo(np.int64).min and c_max < np.iinfo(np.int64).max:
                    df[col] = df[col].astype(np.int64)  
            else:
                if c_min > np.finfo(np.float16).min and c_max < np.finfo(np.float16).max:
                    df[col] = df[col].astype(np.float16)
                elif c_min > np.finfo(np.float32).min and c_max < np.finfo(np.float32).max:
                    df[col] = df[col].astype(np.float32)
                else:
                    df[col] = df[col].astype(np.float64)    
    end_mem = df.memory_usage().sum() / 1024**2
    if verbose: print('Mem. usage decreased to {:5.2f} Mb ({:.1f}% reduction)'.format(end_mem, 100 * (start_mem - end_mem) / start_mem))
    return df

どうやらpandas の各データは静的型付けでいうところの int や long みたいな属性を持っているようで、それを必要十分な容量を持つ型に変換することでメモリの使用量を抑える関数のようです。 for col in df.columns で各行を見て、各行の型、最大値最小値に合わせて適切な型に変換してやる、ということらしい.

## REducing memory
train_df = reduce_mem_usage(train_df)
test_df = reduce_mem_usage(test_df)

weather_train_df = reduce_mem_usage(weather_train_df)
weather_test_df = reduce_mem_usage(weather_test_df)
building_meta_df = reduce_mem_usage(building_meta_df)

この部分のログを見ればわかるように、50%くらいはメモリの使用量が抑えられるらしい。かなりやるべき操作のような気がしますね。

train_df.head()

train_df の頭の方のデータをパッと表示するようです。とりあえず見てみるって感じなのかな

train_df.columns.values

df の行一覧、さらにその値の一覧のarray を計算する感じ

for key, d in train_df.groupby('meter_reading'):
    break
    d.head()
plt.figure(figsize = (20,5))
d['meter'].plot()

正直この部分はよくわからない。多分最後に見た d の["meter"]をプロットする処理になっているはずだが、多分そういう処理ではないと思うし、おそらくd["meter"].plot()が中に入ってくるんじゃないかなと思う

plt.figure(figsize=(8,6))
plt.scatter(range(train_df.shape[0]), np.sort(train_df['meter_reading'].values))
plt.xlabel('index', fontsize=12)
plt.ylabel('meter_reading', fontsize=12)
plt.title("Target Distribution", fontsize=14)
plt.show()

figsize=(8,6)でグラフのサイズを指定して、scatter で散布図をプロットしているみたいですね。ソートしてからプロットしているのでそれは散布図である必要があるのかな？という気はします

plt.figure(figsize = (15,5))
train_df['meter_reading'].plot()

'meter_reading' の値をプロットしてるみたいですね。
plt.figure() でプロットの内容をリセットできるのかな？横軸は指定しない風なので多分indexが自動的に横軸になるのかなあと思います

train_df['meter_reading'].plot(kind='hist',
                            bins=25,
                            figsize=(15, 5),
                           title='Distribution of Target Variable (meter_reading)')
plt.show()

plot()でkind="hist"　と指定するとヒストグラムが作れるらしくて、特にbins=25とか指定すると横に何個分割するかを指定できるようです。例えば階級値が1~1000まであって、bins=25 と指定してやると、40このクラスに分けてカウントするみたいな感じだと思っています

# Load data
train = train_df.set_index(['timestamp'])

# Plot missing values per building/meter
f,a=plt.subplots(1,4,figsize=(20,30))
for meter in np.arange(4):
    df = train[train.meter==meter].copy().reset_index()
    df['timestamp'] = pd.to_timedelta(df.timestamp).dt.total_seconds() / 3600
    df['timestamp'] = df.timestamp.astype(int)
    df.timestamp -= df.timestamp.min()
    missmap = np.empty((1449, df.timestamp.max()+1))
    missmap.fill(np.nan)
    for l in df.values:
        if l[2]!=meter:continue
        missmap[int(l[1]), int(l[0])] = 0 if l[3]==0 else 1
    a[meter].set_title(f'meter {meter:d}')
    sns.heatmap(missmap, cmap='Paired', ax=a[meter], cbar=False)

pd.to_timedelta がよく話からなったのですが、色々な文字列での表現での時間を、一律でtimedelta というオブジェクトに変換して (いい感じに整形して) 返してくれる関数みたいです。それでパースしたものの.dt.total_sconds()がよくわからないですけど、多分秒で表したら何秒か？っていうのを取得しているんだと思います. 60*60で割って時間に直しているみたいなので(total_hours() みたいなのはないのかな？)

missmap = np.empty((1449, df.timestamp.max()+1))
    missmap.fill(np.nan)

この辺でnan で埋めて、

    for l in df.values:
        if l[2]!=meter:continue
        missmap[int(l[1]), int(l[0])] = 0 if l[3]==0 else 1

この辺でbuilding_id, time_stamp を見ながら値が入ってたら1、そうでなければ0 とかで埋めてるって感じ

total = train_df.isnull().sum().sort_values(ascending = False)
percent = (train_df.isnull().sum()/train_df.isnull().count()*100).sort_values(ascending = False)
missing__train_data  = pd.concat([total, percent], axis=1, keys=['Total', 'Percent'])
missing__train_data.head(4)

欠損しているデータの割合がどんな感じなのかを調べる。 .count()でTrue, False 関係なく数えるっていう感じなのだと思います

# Find correlations with the target and sort
correlations = train_df.corr()['meter_reading'].sort_values()

# Display correlations
print('Most Positive Correlations:\n', correlations.tail(15))
print('\nMost Negative Correlations:\n', correlations.head(15))

train_df.corr()で各カラム間の相関係数を計算してくれるようです. そのうちの'meter_reading'の部分をソートして表示しています．

corrs = train_df.corr()


plt.figure(figsize = (20, 8))

# Heatmap of correlations
sns.heatmap(corrs, cmap = plt.cm.RdYlBu_r, vmin = -0.25, annot = True, vmax = 0.6)
plt.title('Correlation Heatmap');

相関係数を計算した後、ヒートマップを表示しています。各引数の意味は、

cmap -> 値に対する色のマッピングを指定.
vmin, vmax -> ヒートマップの最小値、最大値を指定（使い道がよくわからない)
annot -> ヒートマップ上に数値を表示する

だそうです. 色々あって覚えるのが大変そうだなぁー

def plot_dist_col(column):
    '''plot dist curves for train and test weather data for the given column name'''
    fig, ax = plt.subplots(figsize=(10, 10))
    sns.distplot(weather_train_df[column].dropna(), color='green', ax=ax).set_title(column, fontsize=16)
    sns.distplot(weather_test_df[column].dropna(), color='purple', ax=ax).set_title(column, fontsize=16)
    plt.xlabel(column, fontsize=15)
    plt.legend(['train', 'test'])
    plt.show()

かなり重要そうだなという感じがする関数. sns.distplot で分布をプロットしてくれるようです。タイトル付きでプロットした後、プロットにタイトルつけて表示っていう感じかな

ts=train_df.groupby(["timestamp"])["meter_reading"].sum()
ts.astype('float')
plt.figure(figsize=(16,8))
plt.title('meter_reading')
plt.xlabel('timestamp')
plt.ylabel('meter_reading')
plt.plot(ts);

groupby であるカラムを見て同じ値を持つデータに対してまとめて処理を行うっていう感じらしいので、同じタイムスタンプに対して（多分色々なmeter が設定されていて)そのmeter_reading の合計をとるっていう処理になっている？

plt.figure(figsize=(16,6))
plt.plot(ts.rolling(window=12,center=False).mean(),label='Rolling Mean');
plt.plot(ts.rolling(window=12,center=False).std(),label='Rolling sd');
plt.legend();

rolling で、指定した値の数分だけお連続したデータを見て、それに対してmeanやらstdやらで統計量を計算する、という処理らしい。 center=True にすることで（デフォルトはFalseのようだが...) 、窓の真ん中に計算結果の値を置くようになるらしい（デフォルトでは最後に置く). 細かい挙動はよくわかっていないが、なんとなく使い道は想像できる.

import statsmodels.api as sm
# multiplicative
res = sm.tsa.seasonal_decompose(ts.values,freq=12,model="multiplicative")
fig = res.plot()

# Additive model
res = sm.tsa.seasonal_decompose(ts.values,freq=12,model="additive")
fig = res.plot()

statsmodels は統計モデルをパパッと試せるライブラリらしい？そして、seasonal_decomposeは、ある種の時系列データを、周期的なデータ（一週間周期とか)、全体的な傾向を表すデータ(直線的に上がるとか、下がるとかのイメージ？)、そしてランダムなデータに分解してくれるらしい。なんじゃその便利すぎる関数は！（ってかそんなことできるんやって感じだけど...) Additive とか multiplicative とかっていうのは、データを分割するときにその三つの積で表されると思ってやるか、和で表せると思ってやるか、ということの違いのよう。というか、これで説明できるとしたらかなり便利だよなーと思ったり(複雑な変化だと多分自分の目で見てもわからないと思う)

y_mean_time.rolling(window=10).std().plot(figsize=(20, 8))
ax = plt.axhline(y=0.009, color='red')

ローリングしながらstdを計算。plt.axhlineであるyの値で直線を引けるみたいです

```` train_df['meter'] = pd.Categorical(train_df['meter']).rename_categories({0: 'electricity', 1: 'chilledwater', 2: 'steam', 3: 'hotwater'}) daily_train = train_df.copy() daily_train['date'] = daily_train['timestamp'].dt.date daily_train = daily_train.groupby(['date', 'building_id', 'meter']).sum() daily_train

rename_categories で、カテゴリー変数？を変換しているみたいです (0 -> electricity に、みたいな感じ)
dt.date で日付までで日時の情報を丸めて、それを元に groupbyで指定した三つの情報が同じカラムについて合計を取っているよう。なるほどなーという感じですよねー

daily_train_agg = daily_train.groupby(['date', 'meter']).agg(['sum', 'mean', 'idxmax', 'max']) daily_train_agg = daily_train_agg.reset_index() level_0 = daily_train_agg.columns.droplevel(0) level_1 = daily_train_agg.columns.droplevel(1) level_0 = ['' if x == '' else '-' + x for x in level_0] daily_train_agg.columns = level_1 + level_0 daily_train_agg.rename_axis(None, axis=1) daily_train_agg.head()

groupby で date, meter が同じものをまとめて、それに対して sum, mean, idxmax, max を適用.  idx は範囲の中で最大の値を持つような列の名前を返すような関数. 
・・・なんですがちょっとこの部分の処理については理解しきれてないのでおいおい勉強していこうと思っています

fig_total = px.line(daily_train_agg, x='date', y='meter_reading-sum', color='meter', render_mode='svg') fig_total.update_layout(title='Total kWh per energy aspect') fig_total.show()

px.line で折れ線グラフを表示。daily_train_agg を元データとして、'meter'をカラーとして指定するのでそれがキーとして、"meter_reading-sum"の値についてグラフを描画、ということか.

で、どうやらある一つの建物が明らかにおかしなデータを出していることがグラフを見ていると分かるらしい（ほんまか？）。のでそれがどれかを調べていく。

daily_train_agg['building_id_max'] = [x[1] for x in daily_train_agg['meter_reading-idxmax']] daily_train_agg.head()

まず、各日、各メーターについて一番高くなっているidを計算しておく.

def show_building(building, energy_aspects=None): ...

一応定義はされてるけど使われてなかったのでパス

daily_train_electricity = daily_train_agg[daily_train_agg['meter']=='electricity'].copy() daily_train_electricity['building_id_max'] = pd.Categorical(daily_train_electricity['building_id_max']) fig_daily_electricity = px.scatter(daily_train_electricity, x='date', y='meter_reading-max', color='building_id_max', render_mode='svg') fig_daily_electricity.update_layout(title='Maximum consumption values for the day and energy aspect') fig_daily_electricity.show()

` daily_train_agg[daily_train_agg['meter']=='electricity'] ` で ['meter'] == 'electricity' のカラムのidを取ってきて、それを元にコピー.
pd.Categorical でカラムのデータ型をカテゴリー型に変換してくれるらしい. こっちの方がカテゴリーとして扱うときに色々便利なよう.
その次にox.scatter で散布図を作る. x, yで軸にするカラムを指定して, color で色分けに使うカラムを指定. 
可視化された結果を見ると、明らかに値がおかしなことになってるbuilding がいくつかあることがわかりますね．

このエラーが起きている理由は色々と考えられるようですが、何れにせよなんらかの対処をしなければいけない、ということのよう.

ここからが特徴エンジニアリングと、モデリングの部分らしい.

train_df['timestamp'] = pd.to_datetime(train_df['timestamp'])
test_df['timestamp'] = pd.to_datetime(test_df['timestamp'])
weather_train_df['timestamp'] = pd.to_datetime(weather_train_df['timestamp'])
weather_test_df['timestamp'] = pd.to_datetime(weather_test_df['timestamp'])
    
building_meta_df['primary_use'] = building_meta_df['primary_use'].astype('category')

timestamp をto_datetime　でパース.

temp_df = train_df[['building_id']]
temp_df = temp_df.merge(building_meta_df, on=['building_id'], how='left')
del temp_df['building_id']
train_df = pd.concat([train_df, temp_df], axis=1)

temp_df = test_df[['building_id']]
temp_df = temp_df.merge(building_meta_df, on=['building_id'], how='left')

del temp_df['building_id']
test_df = pd.concat([test_df, temp_df], axis=1)
del temp_df, building_meta_df

temp_df = train_df[['site_id','timestamp']]
temp_df = temp_df.merge(weather_train_df, on=['site_id','timestamp'], how='left')

del temp_df['site_id'], temp_df['timestamp']
train_df = pd.concat([train_df, temp_df], axis=1)

temp_df = test_df[['site_id','timestamp']]
temp_df = temp_df.merge(weather_test_df, on=['site_id','timestamp'], how='left')

del temp_df['site_id'], temp_df['timestamp']
test_df = pd.concat([test_df, temp_df], axis=1)

del temp_df, weather_train_df, weather_test_df

building_idやらsite_id, timestamp でDataFrame をマージし、共通する部分を削除してから繋げることで、meta_df を繋げたDFを計算しています

train_df.to_pickle('train_df.pkl')
test_df.to_pickle('test_df.pkl')
   
del train_df, test_df
gc.collect()

pickle によってデータを一時ファイルとして保存することができている

train_df['age'] = train_df['year_built'].max() - train_df['year_built'] + 1
test_df['age'] = test_df['year_built'].max() - test_df['year_built'] + 1

建てられてから何年経っているか？をカラムとして追加する

le = LabelEncoder()
# train_df['primary_use'] = train_df['primary_use'].astype(str)
train_df['primary_use'] = le.fit_transform(train_df['primary_use']).astype(np.int8)

# test_df['primary_use'] = test_df['primary_use'].astype(str)
test_df['primary_use'] = le.fit_transform(test_df['primary_use']).astype(np.int8)

primary_use のカラムをラベル変数に変換する。各ラベルに順序がある場合はその順序を考慮したラベル付け(大きい順とか) にした方が良いらしい.

ここからは、欠損値の扱いについて~~ 欠損値がなぜ欠損しているかを考えることはとても重要であるらしい. 　その埋め方はカーネルに色々書いてありました. 今回のカーネルでは、たくさん値がなくなっているような4つの特徴について、-999で埋めるという手法を取っています

train_df['floor_count'] = train_df['floor_count'].fillna(-999).astype(np.int16)
test_df['floor_count'] = test_df['floor_count'].fillna(-999).astype(np.int16)

train_df['year_built'] = train_df['year_built'].fillna(-999).astype(np.int16)
test_df['year_built'] = test_df['year_built'].fillna(-999).astype(np.int16)

train_df['age'] = train_df['age'].fillna(-999).astype(np.int16)
test_df['age'] = test_df['age'].fillna(-999).astype(np.int16)

train_df['cloud_coverage'] = train_df['cloud_coverage'].fillna(-999).astype(np.int16)
test_df['cloud_coverage'] = test_df['cloud_coverage'].fillna(-999).astype(np.int16)

その次に、たくさんtimestamp に関わる特徴量を計算しています.

drop_cols = ["precip_depth_1_hr", "sea_level_pressure", "wind_direction", "wind_speed","timestamp"]
target = np.log1p(train_df["meter_reading"])  
del train["meter_reading"]
train_df = train_df.drop(drop_cols, axis=1)

ログスケールに直したり、カラムを削除したりしているようです.

あとはlgbを使って学習を回して、っていう感じですね！

感想

とりあえず一通り読んではみましたが、まだまだ自分でかけるという感じではない気がするのでもっと精進していかないといけないですね！

Island Babies

いろんなことを書きます

EDA Kernel を読んでみる

はじめに

感想