◇◆機械学習モデリング前編 Jane Street Market Predictionで市場予測をしよう！◆◇

モチベーション

本職でずっと投資の研究をしているのですが、アイデアや戦略は結構思いついたりします。しかし、経験者であれば分かる通り、上がるか下がるかを予測することは非常に難易度が高いです。

また、アイデア自体は思いついても正直枯渇します。スポンジのように知識を吸収しては、絞り切ってを繰り返すと思います。

なので、ある程度、市場予測の一連の作業を楽に出来ないか考えた所、学習コストはかかるものの、モデリングを鍛えるという結論に到達しました。

本記事では、Kaggleのnotebookを題材にモデリングについて解説していきたいと思います。

ただ、もちろん時系列データの扱いは非常に繊細なので、検証段階で将来の情報リーク等によって良い結果を出すこと等もあり、そういった時系列データの取り扱いの細かい部分は ◆アルゴリズム取引　～投資の定量分析の犯しやすいミスとその解決策～◆ で解説してますので、そちらをご参照くださると助かります。

方針

KaggleのJane Street Market Predictionからnotebookを参考にしながら分析手法やモデリング手法について学んでいきたいと思います。

本記事で採用するnotebookはJane_street_Extensive_EDAで、こちらでは本コンペのデータの理解を深めるためにデータを細かく見ます。

ですので、モデリング編では、データ探索をして、データへの理解度を深めてから、様々な機械学習モデルを用いて市場予測していく方針です。

データの説明

参考:https://www.kaggle.com/c/jane-street-market-prediction/data

このデータセットには、実際の株式市場データを表す匿名化された一連の特徴量 feature_ {0 ... 129} が含まれています。データセットの各行は取引機会を表しており、action を予測します。1は取引を行い、0は取引を渡します。各取引には関連するweightとrespがあり、これらを合わせて取引の収益を表します。date 列は取引の日を表す整数であり、ts_id は時間の順序を表します。匿名化された特徴量に加えて、features.csvでは特徴量に関するメタデータが提供されます。

訓練データtrain.csvには、resp値と、さまざまな期間のリターンを表す他のいくつかのresp_ {1,2,3,4}値が提供されます。これらの変数はテストセットに含まれていません。weight = 0の取引は、完全を期すために意図的にデータセットに含まれていますが、そのような取引はスコア評価に貢献しません。

これは、時系列APIに依存して、モデルが時間的に前向きにならないようにするコード競争です。 APIを使用するには、評価ページの指示に従ってください。ノートブックを送信すると、目に見えないテストで再実行されます。

• コンペのモデルトレーニングフェーズ中、この目に見えないテストセットは、約100万行の履歴データで構成されます。
• ライブ予測フェーズでは、テストセットは定期的に更新されるライブ市場データを使用します。

コンペの第2（予測）フェーズでは、ノートブックの時間制限は、テストセットで提示された取引の数に比例することに注意してください。詳細については、コード要件を参照してください。

成績評価基準

この大会は、効用スコア(utility score)で評価されます。テストセットの各行は、アクション値を予測する取引機会を表します。1は取引を行い、0は取引を見送ります。各トレードjには、リターンを表す重み(weight)と各アセットのリターン(resp)が関連付けられています。

各dateに対して、

\displaystyle p_i = \sum_{j}(weight_{ij} * resp_{ij} * action_{ij})

\displaystyle t = \frac{\sum p_i}{\sqrt{\sum p_i^2}} * \sqrt{\frac{250}{|i|}}

ただし、iはテスト期間における独自な日付の個数です。効用スコアは次のように定義します。

u = \min( \max(t, 0), 6)\sum p_i

解説: ある銘柄jを考えますと株価からリターンを計算出来ます。ここで、株価のリターンをそれぞれ$resp_1, resp_2, ..., resp_i$となっています。また、重みは本コンペで定義された重みですが、ポートフォリオ理論的には銘柄iに資産をどの程度分配したかを表します。なので、本コンペがweight = 0と書いてある銘柄に関してはそもそも評価に貢献しません。

なので、weight * respで個別銘柄の戦略による損益で、actionはそれを取引するかどうかです。$p_i$では全銘柄の投資戦略の損益を表していると言えます。また、$t$に関しては、投資戦略のシャープレシオを表していまして、$\sqrt(250/i)$は投資戦略によって選んだ銘柄の個数を調整するための掛け目です。

最後に$u$ですが、こちらが評価関数になります。意味はシャープレシオ*リターンで、シャープレシオには、上限値と下限値を設けております。

前置きが長くなりましたが、さっそくデータについて見ていきましょう！

notebook解説

インポート

import os
import gc
import sys

import numpy as np
import pandas as pd

import matplotlib.pyplot as plt
import matplotlib.gridspec as gridspec
import seaborn as sns
import plotly.express as px
import plotly.graph_objects as go
%matplotlib inline
sns.set_style('darkgrid')
pd.set_option('display.max_columns', 140)
# from sklearn.preprocessing import StandardScaler as scale
# from sklearn.decomposition import PCA
# from sklearn.neural_network import MLPClassifier as mlp
# from sklearn.linear_model import LogisticRegression as log
# from sklearn.model_selection import train_test_split as split 
# from sklearn.model_selection import GridSearchCV as Grid
# from sklearn.naive_bayes import GaussianNB as GNB
# from sklearn.svm import SVC
# from sklearn.linear_model import LinearRegression as linear
# import xgboost as xgb
# from sklearn.decomposition import IncrementalPCA as ipca
# from sklearn.metrics import (roc_auc_score, precision_score, recall_score, f1_score,
#                              confusion_matrix, accuracy_score, roc_curve, auc)

コメントアウトしている行に関しましては、基本的にモデリングに使う行ですが、今回はデータセットの解説から入るため、こちらは一旦不要です。

データの読み込み

df = pd.read_csv('../input/jane-street-market-prediction/train.csv')

# 元データのサイズを保存
org_len = len(df)

データの大きさを確認

df.info()

<class 'pandas.core.frame.DataFrame'>
RangeIndex: 2390491 entries, 0 to 2390490
Columns: 138 entries, date to ts_id
dtypes: float64(135), int64(3)
memory usage: 2.5 GB

今回の訓練用のデータセットは138列、2,390,491行あります。また、メモリは2.5GB使うため非常に重いです。このままデータを用いるとコンピュータに負担がかかるのでメモリを節約しましょう。

まずはいくつかの変数の型を変更します。

df.feature_0 = df.feature_0.astype(np.int8)
df.date= df.date.astype(np.int16)
df.ts_id = df.ts_id.astype(np.int32)

次のコードは、データフレームをループしながら、-0.0001から0.0001の間のデータの量が非常に少ない場合（0.1％の非常に控えめなスレッショルド）にのみfloat64列をfloat32に変更して、小さな値の列の精度を損なうことを回避します。

for i in df:
    if df[i].dtype == np.float64:
        if (((df[i] < .0001) & (df[i] > -.0001)).mean()) > .001:
            print(i)

weight
resp_1
resp_2
resp_3
resp_4
resp

ここで、これらのカラムに対してメモリ削減のために、型を変更、および、ガベージコレクション(gc)でメモリ管理をします。

for i in df:
    if df[i].dtype == np.float64:
        if (((df[i] < .0001) & (df[i] > -.0001)).mean()) < .001:
            df[i] = df[i].astype(np.float32)
            gc.collect();

では、もう一度、このデータセットのメモリを見てみましょう。

df.info()

<class 'pandas.core.frame.DataFrame'>
RangeIndex: 2390491 entries, 0 to 2390490
Columns: 138 entries, date to ts_id
dtypes: float32(129), float64(6), int16(1), int32(1), int8(1)
memory usage: 1.3 GB

なんと、2.5GBのメモリを消費するデータセットを1.3GBまでに落とすことができました。これで計算も速くなりますし、メモリ不足の問題も多少は改善されると思います。

次に、データセットは時系列データなので、時系列にちゃんとなっているように並び替えましょう。

df.sort_values(by= ['date','ts_id'],inplace=True)

次にactionを定義します。これは、銘柄を買うか売るかでしたね。

df['action'] = np.where(df['resp'] > 0,1,0)
df.action = df.action.astype('category')

ここでnp.whereが出てきましたが、これはnp.where(condition, x, y)で第一引数に条件文を入れて、それが真なら第二引数の値、偽であれば第三引数の値になります。

これはカテゴリカル変数として定義しておきます。つまり、1か0には数値的な意味はなく、買ったか買っていないかを分類するための数字でしかないということです。

EDA(探索的データ解析)

EDAはデータセットを調査して、データの中身について理解を深めることです。これによって、データを用いて予測する際に、どの変数を用いるか、もしくはどのような変数を生成すれば良いのか等について知見を深めることが出来ます。名前はかっこいいですね。枝と読みたい方もいらっしゃるでしょうけど、残念ながら巷ではイーディーエーと読みます。

一般的なデータ分析

Respデータ分析

fig = plt.figure(figsize=(16,6))
ax = plt.subplot(1,1,1)
df.groupby('date')[['resp_1', 'resp_2', 'resp_3', 'resp_4', 'resp']].sum().cumsum().plot(ax=ax)
plt.title('Cumulative Sum of Different RESP\'s',fontsize=18)
plt.xlabel('Date',fontsize=14)
plt.axvspan(0,92,linestyle=':',linewidth=2,label='first 92 days',color='darkorange',alpha=.2)
plt.legend(fontsize=12,ncol=3,loc=2);

最初の92日間でより多くの利益があったことに気付くことができます。これは、多く分析者に、この時点より前に観測を取り除くことをお勧めするという洞察につながります。
また、resp_4の累積合計が最大であるのに対し、resp_1の累積合計は最小であることがわかります。

では、Respの日毎の平均をプロットしてみましょう。

fig = px.line(df.groupby('date')[['resp_1', 'resp_2', 'resp_3', 'resp_4','resp']].mean(),
              x= df.groupby('date')[['resp_1', 'resp_2', 'resp_3', 'resp_4','resp']].mean().index,
              y= ['resp_1', 'resp_2', 'resp_3', 'resp_4','resp'],
              title= '\naverage Resp per day')
fig.layout.xaxis.title = 'Day' 
fig.layout.yaxis.title = 'Avg Resp'
fig.show()

fig,((ax11,ax12,ax13),(ax21,ax22,ax23),(ax31,ax32,ax33),(ax41,ax42,ax43),(ax51,ax52,ax53)) = plt.subplots(5,3,figsize=(18,24))
plt.subplots_adjust(hspace=0.35)
ax11.hist(df.resp,bins=150, color='darkblue',alpha=.6)
ax11.axvline(df.resp.mean()+df.resp.std(),color='darkorange',linestyle=':',linewidth=2)
ax11.axvline(df.resp.mean()-df.resp.std(),color='darkorange',linestyle=':',linewidth=2)
df.resp.plot.hist(bins= 150,ax=ax12,color='darkblue',alpha=.6)
ax12.axvline(df.resp.mean()+df.resp.std(),color='darkorange',linestyle=':',linewidth=2)
ax12.axvline(df.resp.mean()-df.resp.std(),color='darkorange',linestyle=':',linewidth=2)
ax12.set_xlim(-.08,.08)
ax13.hist(df.resp,bins=150, color='darkblue',alpha=.6)
ax13.set_yscale('log')
skew= round(df.resp.skew(),4)
kurt= round(df.resp.kurtosis())
std1= round((((df.resp.mean()-df.resp.std()) < df.resp ) & (df.resp < (df.resp.mean()+df.resp.std()))).mean()*100,2)
props = dict(boxstyle='round', facecolor='white', alpha=0.5)
ax11.text(.02,.96,'μ = {}\nstd = {}\nskewness = {}\nkurtosis = {}\n% values in 1 std = {}%'.format(round(df.resp.mean(),4),round(df.resp.std(),4),skew,kurt,std1),
         transform=ax11.transAxes, verticalalignment='top',bbox=props,fontsize=10)
ax11.set_title('Resp Hist Normal scale',fontsize=14)
ax12.set_title('Resp Hist normal scale zoomed',fontsize=14)
ax13.set_title('Resp Hist with freq on a log scale',fontsize=14);
ax11.set_xlabel('')
ax11.set_ylabel('')
ax12.set_xlabel('')
ax12.set_ylabel('')
ax13.set_xlabel('')
ax13.set_ylabel('')
ax21.hist(df.resp_1,bins=150,color='darkblue',alpha=.6)
ax21.axvline(df.resp_1.mean()+df.resp_1.std(),color='darkorange',linestyle=':',linewidth=2)
ax21.axvline(df.resp_1.mean()-df.resp_1.std(),color='darkorange',linestyle=':',linewidth=2)
df.resp_1.plot.hist(bins= 150,ax=ax22,color='darkblue',alpha=.6)
ax22.axvline(df.resp_1.mean()+df.resp_1.std(),color='darkorange',linestyle=':',linewidth=2)
ax22.axvline(df.resp_1.mean()-df.resp_1.std(),color='darkorange',linestyle=':',linewidth=2)
ax22.set_xlim(-.08,.08)
ax23.hist(df.resp_1,bins=150,color='darkblue',alpha=.6)
ax23.set_yscale('log')
skew= round(df.resp_1.skew(),4)
kurt= round(df.resp_1.kurtosis())
std1= round((((df.resp_1.mean()-df.resp_1.std()) < df.resp_1 ) & (df.resp_1 < (df.resp_1.mean()+df.resp_1.std()))).mean()*100,2)
ax21.text(.02,.96,'μ = {}\nstd = {}\nskewness = {}\nkurtosis = {}\n% values in 1 std = {}%'.format(round(df.resp_1.mean(),4),round(df.resp_1.std(),4),skew,kurt,std1),
         transform=ax21.transAxes, verticalalignment='top',bbox=props,fontsize=10)
ax21.set_title('Resp_1 Hist Normal scale',fontsize=14)
ax22.set_title('Resp_1 Hist normal scale zoomed',fontsize=14)
ax23.set_title('Resp_1 Hist with freq on a log scale',fontsize=14);
ax21.set_xlabel('')
ax21.set_ylabel('')
ax22.set_xlabel('')
ax22.set_ylabel('')
ax23.set_xlabel('')
ax23.set_ylabel('')
ax31.hist(df.resp_2,bins=150,color='darkblue',alpha=.6)
ax31.axvline(df.resp_2.mean()+df.resp_2.std(),color='darkorange',linestyle=':',linewidth=2)
ax31.axvline(df.resp_2.mean()-df.resp_2.std(),color='darkorange',linestyle=':',linewidth=2)
df.resp_2.plot.hist(bins= 150,ax=ax32,color='darkblue',alpha=.6)
ax32.axvline(df.resp_2.mean()+df.resp_2.std(),color='darkorange',linestyle=':',linewidth=2)
ax32.axvline(df.resp_2.mean()-df.resp_2.std(),color='darkorange',linestyle=':',linewidth=2)
ax32.set_xlim(-.08,.08)
ax33.hist(df.resp_2,bins=150, color='darkblue',alpha=.6)
ax33.set_yscale('log')
skew= round(df.resp_2.skew(),4)
kurt= round(df.resp_2.kurtosis())
std1= round((((df.resp_2.mean()-df.resp_2.std()) < df.resp_2 ) & (df.resp_2 < (df.resp_2.mean()+df.resp_2.std()))).mean()*100,2)
ax31.text(.02,.96,'μ = {}\nstd = {}\nskewness = {}\nkurtosis = {}\n% values in 1 std = {}%'.format(round(df.resp_2.mean(),4),round(df.resp_2.std(),4),skew,kurt,std1),
         transform=ax31.transAxes, verticalalignment='top',bbox=props,fontsize=10)
ax31.set_title('Resp_2 Hist Normal scale',fontsize=14)
ax32.set_title('Resp_2 Hist normal scale zoomed',fontsize=14)
ax33.set_title('Resp_2 Hist with freq on a log scale',fontsize=14);
ax31.set_xlabel('')
ax31.set_ylabel('')
ax32.set_xlabel('')
ax32.set_ylabel('')
ax33.set_xlabel('')
ax33.set_ylabel('')
ax41.hist(df.resp_3, color='darkblue',alpha=.6,bins=150)
ax41.axvline(df.resp_3.mean()+df.resp_3.std(),color='darkorange',linestyle=':',linewidth=2)
ax41.axvline(df.resp_3.mean()-df.resp_3.std(),color='darkorange',linestyle=':',linewidth=2)
df.resp_3.plot.hist(bins=150, color='darkblue',alpha=.6,ax=ax42)
ax42.axvline(df.resp_3.mean()+df.resp_3.std(),color='darkorange',linestyle=':',linewidth=2)
ax42.axvline(df.resp_3.mean()-df.resp_3.std(),color='darkorange',linestyle=':',linewidth=2)
ax42.set_xlim(-.08,.08)
ax43.hist(df.resp_3, color='darkblue',alpha=.6,bins=150)
ax43.set_yscale('log')
skew= round(df.resp_3.skew(),4)
kurt= round(df.resp_3.kurtosis())
std1= round((((df.resp_3.mean()-df.resp_3.std()) < df.resp_3 ) & (df.resp_3 < (df.resp_3.mean()+df.resp_3.std()))).mean()*100,2)
ax41.text(.02,.96,'μ = {}\nstd = {}\nskewness = {}\nkurtosis = {}\n% values in 1 std = {}%'.format(round(df.resp_3.mean(),4),round(df.resp_3.std(),4),skew,kurt,std1),
         transform=ax41.transAxes, verticalalignment='top',bbox=props,fontsize=10)
ax41.set_title('Resp_3 Hist Normal scale',fontsize=14)
ax42.set_title('Resp_3 Hist normal scale zoomed',fontsize=14)
ax43.set_title('Resp_3 Hist with freq on a log scale',fontsize=14);
ax41.set_xlabel('')
ax41.set_ylabel('')
ax42.set_xlabel('')
ax42.set_ylabel('')
ax43.set_xlabel('')
ax43.set_ylabel('')
ax51.hist(df.resp_4,bins=150, color='darkblue',alpha=.6)
ax51.axvline(df.resp_4.mean()+df.resp_4.std(),color='darkorange',linestyle=':',linewidth=2)
ax51.axvline(df.resp_4.mean()-df.resp_4.std(),color='darkorange',linestyle=':',linewidth=2)
df.resp_4.plot.hist(bins= 150,color='darkblue',alpha=.6,ax=ax52)
ax52.axvline(df.resp_4.mean()+df.resp_4.std(),color='darkorange',linestyle=':',linewidth=2)
ax52.axvline(df.resp_4.mean()-df.resp_4.std(),color='darkorange',linestyle=':',linewidth=2)
ax52.set_xlim(-.08,.08)
ax53.hist(df.resp_4,bins=150,color='darkblue',alpha=.6)
ax53.set_yscale('log')
skew= round(df.resp_4.skew(),4)
kurt= round(df.resp_4.kurtosis())
std1= round((((df.resp_4.mean()-df.resp_4.std()) < df.resp_4 ) & (df.resp_4 < (df.resp_4.mean()+df.resp_4.std()))).mean()*100,2)
ax51.text(.02,.96,'μ = {}\nstd = {}\nskewness = {}\nkurtosis = {}\n% values in 1 std = {}%'.format(round(df.resp_4.mean(),4),round(df.resp_4.std(),4),skew,kurt,std1),
         transform=ax51.transAxes, verticalalignment='top',bbox=props,fontsize=10)
ax51.set_title('Resp_4 Hist Normal scale',fontsize=14)
ax52.set_title('Resp_4 Hist normal scale zoomed',fontsize=14)
ax53.set_title('Resp_4 Hist with freq on a log scale',fontsize=14)
ax51.set_xlabel('')
ax51.set_ylabel('')
ax52.set_xlabel('')
ax52.set_ylabel('')
ax53.set_xlabel('')
ax53.set_ylabel('')
fig.suptitle('RESPs Historgrams on Different Scales',fontsize=18,y=.92);

次にペアプロットです。

sns.pairplot(df[['resp_1', 'resp_2', 'resp_3', 'resp_4', 'resp']],corner=True);

この図に示されている関係と、この前の図に示されている標準偏差と分布から、RespはResp_4との関連性が高く、Resp_1とResp_2は相互に関連性が高いことがわかります。より長い期間の投資はより多くのリターンとより高いリスクに繋がります。

投資期間に関する基礎知識

投資期間は、特定の目標のために投資を行うことが期待される期間です。投資は一般に、株式（リスクが高い）と債券（リスクが低い）の2つの主要なカテゴリに分類されます。投資期間(以下、タイム・ホライズン)が長ければ長いほど、投資家はより積極的またはリスクの高いポートフォリオを構築できます。タイム・ホライズンが短ければ短いほど、投資家はより保守的またはリスクが少なくなります。

日付(date)

今回は日付がちょうど500日分あるため、これは大体２年分のデータであると判断して良いです。理由は、株式市場の営業日が年間大体250日であるからです。$t$の定義式の掛け目も一年間=250営業日と定義しているように思われます。

では、日別にrespがある銘柄数について見てみましょう。株式市場では、もちろん全銘柄に対してリターンは存在しますが、本コンペでは、ある程度投資対象を絞っているため、全銘柄について調査しなくて済みます。

fig = px.area(data_frame= df.groupby('date')[['resp']].count(),title='Number of operation per day')
fig.update_traces( showlegend = False)
fig.layout.xaxis.title = 'Day' 
fig.layout.yaxis.title = 'Number of operations'
fig.show()

日別合計リターン

fig = px.area(data_frame= df.groupby('date')[['resp']].sum(),title='Resp sum of operation per day')
fig.update_traces( showlegend = False)
fig.layout.xaxis.title = 'Day' 
fig.layout.yaxis.title = 'Resp sum'
fig.show()

respの分散が大きいことが伺えます。

それぞれの平均に対して20営業日（1ヶ月=約20営業日を示唆）の移動標準偏差を作成します。

date_df = df.groupby('date')[['resp']].mean()
std20 = []
for i in range(len(date_df)):
    if i <20:
        std20.append(np.nan)
    else:
        moving_std = date_df['resp'][i-20:i].std()
        std20.append(moving_std)
date_df['moving_std'] = std20
date_df.tail(2)

	resp	moving_std
date
498	0.001554	0.001962
499	0.000299	0.001969

fig = px.line(data_frame=date_df,y=['resp','moving_std'],title='Average Resp & 20 day moving standard deviation')
fig.layout.xaxis.title = 'Day' 
fig.layout.yaxis.title = 'Avg Resp'
fig.show()

ここで、毎日の各respの標準偏差を確認します。

fig, (ax1,ax2) = plt.subplots(2,1,figsize=(14,12))
df.groupby('date')[['resp_1', 'resp_2', 'resp_3', 'resp_4']].std().plot(ax=ax1,color=['steelblue','darkorange','red','green'],alpha=.8)
df.groupby('date')[['resp_1', 'resp_2', 'resp_3', 'resp_4']].std().plot.kde(ax=ax2)
fig.suptitle('Resp\'s Std',fontsize=18,y=.96)
ax2.set_xlabel('')
ax1.set_xlabel('')
ax2.set_title('kde of each resp std', fontsize=14)
ax1.set_title('std of Resp\'s for each trading day',fontsize=14);

前に述べたように、respの標準偏差は、長期的な投資とともに増加するようです。
また、多くのKagglerが80日後に何らかの取引モデルの調整が行われた可能性があると述べたため、最初の100日間で偏差がやや高かったことにも気付くことができます。

重み(weight)

fig = plt.figure(figsize=(18,7))
grid =  gridspec.GridSpec(2,3,figure=fig,hspace=.3,wspace=.2)
ax1 = fig.add_subplot(grid[0, 0])
ax2 = fig.add_subplot(grid[0, 1])
ax3 = fig.add_subplot(grid[1, 0])
ax4 = fig.add_subplot(grid[1, 1])
ax5 = fig.add_subplot(grid[:, 2])
sns.boxplot(x = df.weight,width=.5,ax=ax1)
ax2.hist(df.weight, color='#404788ff',alpha=.6, bins= list([-.05] + list(10**np.arange(-2,2.24,.05))))
ax2.set_xscale('symlog')
ax2.set_xlim(-.05,227)
sns.boxplot(x = df.weight[df.weight != 0],width=.5,ax=ax3)
ax1.set_title('Weights including zero weights',fontsize=14)
ax3.set_title('Weights not including zero weights',fontsize=14)
ax2.set_title('Weights including zero weights (log)',fontsize=14)
ax4.set_title('Weights not including zero weights (log)',fontsize=14)
props = dict(boxstyle='round', facecolor='white', alpha=0.4)
ax1.text(.2,.9,'μ = {}    std = {}\nmin = {}    max = {}'.format(round(df.weight.mean(),3),round(df.weight.std(),3),round(df.weight.min(),3),round(df.weight.max(),3)),
         transform=ax1.transAxes, verticalalignment='top',bbox=props,fontsize=12)
ax3.text(.2,.9,'μ = {}        std = {}\nmin = {}    max = {}'.format(round(df.weight[df.weight != 0].mean(),3),round(df.weight[df.weight != 0].std(),3),
                                                              round(df.weight[df.weight != 0].min(),3),round(df.weight[df.weight != 0].max(),3)),
         transform=ax3.transAxes, verticalalignment='top',bbox=props,fontsize=12)
ax4.hist(df.weight[df.weight !=0],color='#404788ff',alpha=.6,bins=10**np.arange(-2.16,2.24,.05))
ax4.set_xscale('log')
ax4.set_xticks((.01,.03,.1,.3,1,3,10,30,100))
ax4.set_xticklabels((.01,.03,.1,.3,1,3,10,30,100))
ax5.pie(((df.weight==0).mean(),(1-(df.weight==0).mean())),startangle=300,wedgeprops=dict(width=0.5),
        labels=('Zeros\n{}%'.format(round((df.weight==0).mean()*100,2)),'Nonzeros\n{}%'.format(round((1-(df.weight==0).mean())*100,2))),
        textprops={'fontsize': 12},colors=['#404788ff','#55c667ff'])
ax5.set_title('Zeros vs non-zero weights',fontsize=14)
ax1.set_xlabel('')
ax2.set_xlabel('')
ax3.set_xlabel('')
ax2.set_ylabel('')
ax5.set_ylabel('')
ax4.set_xlabel('');

fig = plt.figure(figsize=(15,10))
fig.suptitle('Nonzero weights histogram in different scales',fontsize=18)
ax1 = plt.subplot(3,1,1)
ax1.hist(df.weight[df.weight !=0],color='darkblue',alpha=.7, bins=10**np.arange(-2.16,2.23,.05))
plt.xscale('log')
plt.xticks((.01,.03,.1,.3,1,3,10,30,100),(.01,.03,.1,.3,1,3,10,30,100))
ax2 = plt.subplot(3,1,2)
sns.distplot(df.weight[df.weight != 0], color='darkblue', bins=400, ax=ax2) 
ax3 = plt.subplot(3,1,3)
ax3.hist(df.weight[(df.weight !=0) & (df.weight < 3.197 )],color='darkblue',alpha=.7, bins=200)
ax3.set_xlim(0,3.3)
ax2.set_xlabel('') 
ax1.set_title('All values (log-scale)',fontsize=14)
ax2.set_title('kde of the distribution',fontsize=14)
ax3.set_title('75% of the Values',fontsize=14)
plt.subplots_adjust(hspace=.4);

重みの値が0に集中していて、かなり分散が激しいですね。特に右端の外れ値に関してより詳しく見ていきましょう。

fig, (ax1,ax2) = plt.subplots(2,1,figsize=(16,8))
fig.suptitle('Weight outliers',fontsize=18)
sns.boxplot(df.weight,width=.5, ax=ax1)
ax1.axvline(np.percentile(df.weight,95), color= 'green',label='95.0%',linestyle=':',linewidth=3)
ax1.axvline(np.percentile(df.weight,99), color= 'darkblue',label='99.0%',linestyle=':',linewidth=3)
ax1.axvline(np.percentile(df.weight,99.9), color= 'darkorange',label='99.9%',linestyle=':',linewidth=3)
ax1.axvline(np.percentile(df.weight,99.99), color= 'magenta',label='99.99%',linestyle=':',linewidth=3)
ax1.legend(fontsize=13)
sns.boxplot(df.weight[df.weight !=0],width=.5, ax=ax2)
ax2.axvline(np.percentile(df.weight[df.weight !=0],95), color= 'green',label='95.0%',linestyle=':',linewidth=3)
ax2.axvline(np.percentile(df.weight[df.weight !=0],99), color= 'darkblue',label='99.0%',linestyle=':',linewidth=3)
ax2.axvline(np.percentile(df.weight[df.weight !=0],99.9), color= 'darkorange',label='99.9%',linestyle=':',linewidth=3)
ax2.axvline(np.percentile(df.weight[df.weight !=0],99.99), color= 'magenta',label='99.99%',linestyle=':',linewidth=3)
ax2.legend(fontsize=13)
ax1.set_title('All weights', fontsize= 14)
ax2.set_title('Non-zero weights', fontsize= 14)
ax1.set_xlabel('')
ax2.set_xlabel('');

sns.scatterplot(data=df, x='resp',y='weight', color= 'blue', alpha=.3)
plt.title('Resp vs Weight\ncorrelation={}'.format(round(df.weight.corr(df.resp),4)));

重みはRespと線形に相関していないことがわかりますが、重みが大きいほどResp値が低いことにのみ関連していることは明らかです。

特徴量データ分析

特徴量csvを読み込む

df_f =  pd.read_csv('../input/jane-street-market-prediction/features.csv')
df_f.head(5)

	feature	tag_0	tag_1	tag_2	tag_3	tag_4	tag_5	tag_6	tag_7	tag_8	tag_9	tag_10	tag_11	tag_12	tag_13	tag_14	tag_15	tag_16	tag_17	tag_18	tag_19	tag_20	tag_21	tag_22	tag_23	tag_24	tag_25	tag_26	tag_27	tag_28
0	feature_0	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False
1	feature_1	False	False	False	False	False	False	True	True	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False
2	feature_2	False	False	False	False	False	False	True	True	False	True	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False
3	feature_3	False	False	False	False	False	False	True	False	True	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False
4	feature_4	False	False	False	False	False	False	True	False	True	True	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False	False

fig = px.bar(df_f.set_index('feature').T.sum(), title='Number of tags for each feature')
fig.layout.xaxis.tickangle = 300
fig.update_traces( showlegend = False)
fig.layout.xaxis. dtick = 5
fig.layout.xaxis.title = ''
fig.layout.yaxis.title = ''
fig.show()

欠損値

fig = px.bar(x = df.isnull().sum().index,y= df.isnull().sum().values,title= 'Number of Null values')
fig.layout.xaxis.tickangle = 300
fig.layout.xaxis. dtick = 5
fig.layout.yaxis. dtick = 100000
fig.layout.xaxis.title = ''
fig.layout.yaxis.title = ''
fig.layout.xaxis.showgrid = True
fig.show()

欠損値が多いカラムが目立ちますね。データの10%以上が欠損しているカラムを見つけましょう。

nulls = df.isnull().sum()
nulls_list = list(nulls[nulls >(0.1 * len(df))].index)
nulls_list

['feature_7',
 'feature_8',
 'feature_17',
 'feature_18',
 'feature_27',
 'feature_28',
 'feature_72',
 'feature_78',
 'feature_84',
 'feature_90',
 'feature_96',
 'feature_102',
 'feature_108',
 'feature_114']

null値の数にはある種のパターンがあるため、これらの特徴量間の関係を調べます。※参照notebookでは、色が表示されますが、こちらでは、普通のテーブルを使用します。

df[['resp','resp_1','resp_2','resp_3','resp_4','weight']+nulls_list].corr().style.background_gradient(cmap='coolwarm')

	resp	resp_1	resp_2	resp_3	resp_4	weight	feature_7	feature_8	feature_17	feature_18	feature_27	feature_28	feature_72	feature_78	feature_84	feature_90	feature_96	feature_102	feature_108	feature_114
resp	1.000000	0.452159	0.595274	0.815972	0.956197	-0.006948	0.021207	0.014660	-0.028277	-0.017448	0.048293	0.032502	-0.002034	0.000072	0.011500	0.003604	0.013469	0.004613	0.013209	0.004385
resp_1	0.452159	1.000000	0.890214	0.678381	0.358949	-0.016091	0.033919	0.026142	0.004333	0.006632	0.036125	0.025671	-0.000889	0.001371	0.008330	0.001973	0.010467	0.002994	0.010255	0.002969
resp_2	0.595274	0.890214	1.000000	0.823383	0.469230	-0.012955	0.030800	0.025124	-0.005811	0.000202	0.039829	0.028633	-0.000642	0.001266	0.008487	0.001454	0.010247	0.002223	0.010270	0.002225
resp_3	0.815972	0.678381	0.823383	1.000000	0.805952	-0.009229	0.019317	0.014859	-0.002078	0.000566	0.023373	0.016951	-0.004490	-0.003221	0.007704	0.000655	0.011542	0.003766	0.010231	0.002438
resp_4	0.956197	0.358949	0.469230	0.805952	1.000000	-0.005817	0.015059	0.009246	-0.018851	-0.012306	0.033533	0.022307	-0.003965	-0.002240	0.009916	0.003087	0.012989	0.005355	0.011978	0.004442
weight	-0.006948	-0.016091	-0.012955	-0.009229	-0.005817	1.000000	-0.016164	-0.000819	-0.011541	-0.002568	-0.010643	0.001118	0.002058	0.000367	-0.041597	0.052486	-0.044057	0.056261	-0.046157	0.053893
feature_7	0.021207	0.033919	0.030800	0.019317	0.015059	-0.016164	1.000000	0.898832	0.402987	0.374951	0.487964	0.430844	-0.009719	-0.006903	-0.015400	-0.054359	-0.009083	-0.051083	-0.001724	-0.046770
feature_8	0.014660	0.026142	0.025124	0.014859	0.009246	-0.000819	0.898832	1.000000	0.373778	0.425576	0.424587	0.463257	-0.003916	0.001188	-0.002604	-0.038638	0.000781	-0.039802	0.010609	-0.034254
feature_17	-0.028277	0.004333	-0.005811	-0.002078	-0.018851	-0.011541	0.402987	0.373778	1.000000	0.906956	-0.501928	-0.464884	-0.070201	-0.069432	-0.021351	-0.053571	0.039123	0.007819	0.016894	-0.020790
feature_18	-0.017448	0.006632	0.000202	0.000566	-0.012306	-0.002568	0.374951	0.425576	0.906956	1.000000	-0.446295	-0.503707	-0.078615	-0.074634	-0.021682	-0.054821	0.040938	0.006212	0.017579	-0.023075
feature_27	0.048293	0.036125	0.039829	0.023373	0.033533	-0.010643	0.487964	0.424587	-0.501928	-0.446295	1.000000	0.900539	0.059322	0.060716	0.020922	-0.004243	-0.029925	-0.059175	0.000032	-0.028572
feature_28	0.032502	0.025671	0.028633	0.016951	0.022307	0.001118	0.430844	0.463257	-0.464884	-0.503707	0.900539	1.000000	0.074591	0.075029	0.030016	0.005652	-0.027775	-0.055675	0.006896	-0.021918
feature_72	-0.002034	-0.000889	-0.000642	-0.004490	-0.003965	0.002058	-0.009719	-0.003916	-0.070201	-0.078615	0.059322	0.074591	1.000000	0.773391	0.268906	0.209429	-0.208541	-0.144612	0.029765	0.032632
feature_78	0.000072	0.001371	0.001266	-0.003221	-0.002240	0.000367	-0.006903	0.001188	-0.069432	-0.074634	0.060716	0.075029	0.773391	1.000000	0.221750	0.234574	-0.149852	-0.183428	0.035854	0.024444
feature_84	0.011500	0.008330	0.008487	0.007704	0.009916	-0.041597	-0.015400	-0.002604	-0.021351	-0.021682	0.020922	0.030016	0.268906	0.221750	1.000000	0.645486	0.873262	0.534054	0.951272	0.608618
feature_90	0.003604	0.001973	0.001454	0.000655	0.003087	0.052486	-0.054359	-0.038638	-0.053571	-0.054821	-0.004243	0.005652	0.209429	0.234574	0.645486	1.000000	0.547641	0.891866	0.628328	0.971499
feature_96	0.013469	0.010467	0.010247	0.011542	0.012989	-0.044057	-0.009083	0.000781	0.039123	0.040938	-0.029925	-0.027775	-0.208541	-0.149852	0.873262	0.547641	1.000000	0.633677	0.956765	0.609289
feature_102	0.004613	0.002994	0.002223	0.003766	0.005355	0.056261	-0.051083	-0.039802	0.007819	0.006212	-0.059175	-0.055675	-0.144612	-0.183428	0.534054	0.891866	0.633677	1.000000	0.618610	0.971214
feature_108	0.013209	0.010255	0.010270	0.010231	0.011978	-0.046157	-0.001724	0.010609	0.016894	0.017579	0.000032	0.006896	0.029765	0.035854	0.951272	0.628328	0.956765	0.618610	1.000000	0.645604
feature_114	0.004385	0.002969	0.002225	0.002438	0.004442	0.053893	-0.046770	-0.034254	-0.020790	-0.023075	-0.028572	-0.021918	0.032632	0.024444	0.608618	0.971499	0.609289	0.971214	0.645604	1.000000

これらのカラムのnull値の数は膨大であるため（25万を超える）、これらの特徴量とRespおよび重みとの相関関係がないため、null値が10％を超える特徴量は取り除きます。

df.drop(columns=nulls_list,inplace=True)

残りの欠損値については、最初に変動係数に注目します。

(df.iloc[:,7:-2].std() / df.iloc[:,7:-2].mean()).head(5)

feature_0    101.635945
feature_1      6.652154
feature_2      6.905722
feature_3    217.145702
feature_4    421.125070
dtype: float64

平均値がゼロに近いため、変動係数は信頼できないようです

これで、特徴量の分布を俯瞰することができます。

df.iloc[:,7:-2].hist(bins=100,figsize=(20,74),layout=(29,4));

次のコードは、特徴量の分布の包括的で確実な理解を得るために、平均のプロットで水平箱ひげ図のグリッドを作成します。

データが強く中央に偏っているため、0.1％から99.9％のひげを使用して極端な外れ値を取り除いたことに注意してください。

fig = plt.figure(figsize=(20,80))
fig.suptitle('Features Box plot with 0.1% 99.9% whiskers',fontsize=22, y=.89)
grid =  gridspec.GridSpec(29,4,figure=fig,hspace=.5,wspace=.05)
featstr = [i for i in df.columns[7:-2]]
counter = 0
for i in range(29):
    for j in range(4):
        subf = fig.add_subplot(grid[i, j]);
        sns.boxplot(x= df[featstr[counter]],saturation=.5,color= 'blue', ax= subf,width=.5,whis=(.1,99.9));
        subf.axvline(df[featstr[counter]].mean(),color= 'darkorange', label='Mean', linestyle=':',linewidth=3)
        subf.set_xlabel('')
        subf.set_title('{}'.format(featstr[counter]),fontsize=16)
        counter += 1
        gc.collect()
plt.show();

各特徴量の分布に影響を与える外れ値がたくさんあることがわかります。

また、値の大部分は平均を中心としているため、平均を使用してnull値を補完します。

df.fillna(df.mean(axis=0),inplace=True)

特徴量の成長

ここでは、日別の平均のプロットの累積和を見て、特徴量の増加について着目します。

df.groupby('date')[featstr].mean().cumsum().plot(layout=(29,4),subplots=True,figsize=(20,82),xlabel='')
fig = plt.gcf()
fig.text(0.5, 0.19, 'Date',ha='center', fontsize = 24)
fig.suptitle('Cumulative features means per day',fontsize=24,y=.886);

多くの特徴量の日別平均の累積和は直線的に増加しているように見えますが、81、82、83などの一部の特徴量は実際には減少しており、特徴量3のように変動する特徴量もあります。

特徴量同士の相関

相関行列のデータフレームを作成します。

corr = df.iloc[:,7:-2].corr()
corr.style.background_gradient(cmap='coolwarm')

この相関行列は非常に大きい(130✖︎130)ため、ここでは表示を割愛します。

相関行列のヒートマップ

fig = plt.figure(figsize=(18,12))
ax = plt.subplot(1,1,1)
sns.heatmap(corr,ax= ax, cmap='coolwarm');

特徴量間に多重共線性がたくさんあるように見えます。また、特徴量空間にカップルのパターンがあるように見えます。このパターンは、feature_41などの一部のフィーチャで壊れています。

featstr2 = [ i for i in featstr if i not in ['feature_41','feature_64']]
len(featstr)

116

fig = plt.figure(figsize=(22,44))
grid =  gridspec.GridSpec(12,5,figure=fig,hspace=.5,wspace=.2)
counter = 1
for i in range(12):
    for j in range(5):
        if counter == 113:
            break
        subf = fig.add_subplot(grid[i, j]);
        sns.scatterplot(x= df[featstr2[counter]], y = df[featstr2[counter+1]], ax= subf);
        cor = round(df[featstr2[counter]].corr(df[featstr2[counter+1]]) * 100,2)
        subf.set_xlabel('')
        subf.set_ylabel('')
        subf.set_title('{} & {}\nCorrelation = {}%'.format(featstr2[counter],featstr2[counter+1],cor),fontsize=14)
        counter += 2
        gc.collect();  

想定通り、これらは財務関連の特徴量であるため、多くの機能は互いに強く相関しています。

次に、相関性の高い特徴量のグループを調査します。

まずは[feature_19, feature_20, feature_21, feature_22, feature_23, feature_24, feature_25, feature_26, feature_29]から始めます。なぜなら、多重共線生が見えるからです。

plt.figure(figsize=(12,6)) 
sns.heatmap(df[featstr2[15:23]].corr(),center=0,cmap='coolwarm',annot=True,cbar=False);

sns.pairplot(df[featstr2[15:23]],corner=True);

ピアソンの相関係数がこれらの特徴間で非常に高いという事実にもかかわらず、関係は完全に線形ではありませんが、外れ値が散布図の形状に影響を与えることにも注意できます。

では、他のグループについても見てみましょう。

plt.figure(figsize=(12,6)) 
sns.heatmap(df[featstr2[23:31]].corr(),center=0,cmap='coolwarm',annot=True,cbar=False);

sns.pairplot(df[featstr2[23:31]],corner=True);

要点は、他のクラスターとほぼ同じように見えますが、これらの特徴量クラスターの両方が互いに負の相関関係にあることにも言及する価値があります。

plt.figure(figsize=(18,6)) 
sns.heatmap(df[featstr2[15:31]].corr(),center=0,cmap='coolwarm',annot=True,cbar=False);

外れ値

まずは特徴量の平均値について見ましょう。

fig = px.bar(df[featstr].mean(), title='Features mean values')
fig.layout.xaxis.tickangle = 300
fig.update_traces(showlegend = False)
fig.layout.xaxis. dtick = 5
fig.layout.xaxis.title = ''
fig.layout.yaxis.title = ''
fig.show()

次に最大値について見ましょう。

fig = px.bar(df[featstr].max(), title='Features Max Values')
fig.layout.xaxis.tickangle = 300
fig.update_traces(showlegend = False)
fig.layout.xaxis. dtick = 5
fig.layout.xaxis.title = ''
fig.layout.yaxis.title = ''
fig.show()

次は最小値です。

fig = px.bar(df[featstr].min(), title='Features Min Values')
fig.layout.xaxis.tickangle = 300
fig.update_traces(showlegend = False)
fig.layout.xaxis. dtick = 5
fig.layout.xaxis.title = ''
fig.layout.yaxis.title = ''
fig.show()

fig, (ax1,ax2,ax3)= plt.subplots(3,1,figsize=(10,12))
plt.subplots_adjust(hspace=.3)
sns.distplot(df[featstr].max(),ax= ax1 )
sns.distplot(df[featstr].min(),ax= ax2)
sns.distplot(df[featstr].mean(),ax= ax3)
fig.suptitle('distribution of mean max and min for features',fontsize=16)
ax1.set_title('distribution  of features max values',fontsize=14)
ax1.text(.82,.56,'std = {}'.format(round(df[featstr].max().std(),2)),transform=ax1.transAxes, verticalalignment='top',bbox=props,fontsize=12)
ax2.set_title('distribution  of features min values',fontsize=14)
ax2.text(.82,.56,'std = {}'.format(round(df[featstr].min().std(),2)),transform=ax2.transAxes, verticalalignment='top',bbox=props,fontsize=12)
ax3.set_title('distribution  of features mean values',fontsize=14)
ax3.text(.82,.56,'std = {}'.format(round(df[featstr].mean().std(),2)),transform=ax3.transAxes, verticalalignment='top',bbox=props,fontsize=12);

平均値は互いにそれほど異ならないという事実にもかかわらず、最小値と最大値は非常に偏った分布で非常にずれていることがわかります。

外れ値へのより統計的な探索

for i in featstr[1:]:
    print('{}\n0.1%:99.9% are between: {}\nmax: {}\nmin: {}\n75% are under: {}'.format(i,np.percentile(df[i],(.1,99.9)), df[i].max(),df[i].min(),np.percentile(df[i],75)),
         '\n===============================')

feature_1
0.1%:99.9% are between: [-3.1720264  16.70666595]
max: 74.42988586425781
min: -3.1720263957977295
75% are under: 1.5784167051315308 
===============================
feature_2
0.1%:99.9% are between: [-3.09318233 16.21772709]
max: 148.07630920410156
min: -3.093182325363159
75% are under: 1.52639901638031 
===============================
feature_3
0.1%:99.9% are between: [-8.56044363  8.42387932]
max: 25.87166404724121
min: -25.424652099609375
75% are under: 1.0469255447387695 
===============================
feature_4
0.1%:99.9% are between: [-7.86522797  7.90286264]
...

df[(df.feature_56== df.feature_56.max())|(df.feature_57== df.feature_57.max())|(df.feature_58== df.feature_58.max()) | (df.feature_59== df.feature_59.max())]

	date	weight	resp_1	resp_2	resp_3	resp_4	resp	feature_0	feature_1	feature_2	feature_3	feature_4	feature_5	feature_6	feature_9	feature_10	feature_11	feature_12	feature_13	feature_14	feature_15	feature_16	feature_19	feature_20	feature_21	feature_22	feature_23	feature_24	feature_25	feature_26	feature_29	feature_30	feature_31	feature_32	feature_33	feature_34	feature_35	feature_36	feature_37	feature_38	feature_39	feature_40	feature_41	feature_42	feature_43	feature_44	feature_45	feature_46	feature_47	feature_48	feature_49	feature_50	feature_51	feature_52	feature_53	feature_54	feature_55	feature_56	feature_57	feature_58	feature_59	feature_60	feature_61	feature_62	feature_63	feature_64	feature_65	feature_66	feature_67	feature_68	feature_69	feature_70	feature_71	feature_73	feature_74	feature_75	feature_76	feature_77	feature_79	feature_80	feature_81	feature_82	feature_83	feature_85	feature_86	feature_87	feature_88	feature_89	feature_91	feature_92	feature_93	feature_94	feature_95	feature_97	feature_98	feature_99	feature_100	feature_101	feature_103	feature_104	feature_105	feature_106	feature_107	feature_109	feature_110	feature_111	feature_112	feature_113	feature_115	feature_116	feature_117	feature_118	feature_119	feature_120	feature_121	feature_122	feature_123	feature_124	feature_125	feature_126	feature_127	feature_128	feature_129	ts_id	action
235877	38	3.339560	-0.001017	-0.007426	-0.024829	-0.046755	-0.033357	-1	5.399380	5.108447	3.613957	3.282637	1.874917	1.627905	3.537130	2.426174	2.572712	3.605390	-3.484599	-2.825172	2.857596	3.233086	-3.352329	-3.912990	-12.886112	-10.339006	-4.377385	-5.277115	-3.872348	-4.283693	3.506302	4.035473	3.110046	4.478292	3.713556	3.731210	4.706809	4.334348	-3.058521	-3.192860	3.088272	4.539015	-0.136437	3.195502	0.780103	1.167637	-0.327845	0.913224	1.520857	1.316261	1.053037	0.099054	0.855302	1.775148	3.746835	0.985716	388.749817	66.636330	660.351135	800.703613	350.911865	-0.671831	-0.569680	-1.486158	-1.611610	-2.942103	-1.607305	-1.071010	-1.565248	-1.505246	2.665544	-0.095570	1.875405	0.898643	0.415992	1.163417	0.861380	0.671799	0.222680	0.190999	0.681430	0.211753	0.399736	1.027329	0.849244	1.615301	0.705474	-1.423562	-0.041847	-0.107219	-0.006064	-0.310370	-2.582340	-1.338859	0.884709	-0.057046	-1.194013	-1.427599	-0.940190	-0.260405	-0.789533	-1.781693	-2.874468	-0.429607	0.747440	0.427909	-0.787168	-1.908855	-1.880670	-0.350742	-0.874343	-2.260690	-1.953951	6.032958	4.065476	3.634447	1.015100	5.885513	4.415474	8.390346	4.077566	5.459621	2.131100	235877	0
282230	44	10.859674	-0.001306	-0.012651	-0.016530	-0.052089	-0.048993	1	7.414509	9.399608	5.123063	6.146243	4.196547	4.889341	11.098932	11.976443	2.739184	4.887201	6.628507	7.351760	8.330588	11.728739	9.664564	13.641744	5.919429	5.534737	6.929136	9.939451	9.581731	13.009244	-2.449863	-3.702712	-1.050895	-2.082486	-1.886110	-2.196391	-3.542542	-3.657356	5.737075	7.694194	-3.350399	-6.203264	-0.113139	1.447529	-0.254701	7.538714	1.115668	2.226885	3.178381	4.724526	6.178386	1.928744	3.180385	-0.291821	3.230325	2.360054	746.739868	254.094452	816.712097	138.990021	55.821522	-1.349714	-1.129695	-0.094330	-0.107096	-2.041466	-0.597949	-0.441602	-1.366051	-1.314805	16.274912	13.676999	19.838009	-5.815048	-2.006737	-3.691047	-5.506912	0.018003	-2.335785	-0.879096	-2.030294	-2.191662	0.140086	-1.238212	-0.142781	-2.321156	-1.641860	-1.393412	-1.515613	-1.216623	-1.746285	-1.086886	-2.817912	9.008026	2.096552	3.640483	6.876660	-1.182548	3.427998	0.083712	1.071575	3.375455	-2.786946	6.497206	0.910018	1.932831	5.401005	-1.713007	2.378500	-0.578150	0.342619	1.677919	-2.013862	9.998623	10.710449	3.118533	1.845357	6.339649	7.125723	4.243608	2.327895	3.403054	1.688110	282230	0
286937	44	0.219883	0.009060	0.000275	0.008128	0.006852	-0.004134	1	2.464552	3.314075	6.717828	8.427411	5.999205	7.184940	7.984714	8.634279	10.861753	16.884825	16.154589	19.937834	22.595144	31.915413	10.789083	15.700766	14.886054	16.765507	17.146122	24.403795	24.710970	34.675587	-7.076792	-10.310453	-5.595434	-9.575839	-9.691500	-12.738782	-16.241947	-20.111063	7.307025	10.078041	-1.336552	-2.614510	0.127066	0.966415	-0.375578	8.362889	1.284996	-1.139717	-1.475705	-1.193366	-0.914533	-2.279577	-1.075379	-1.819851	8.140661	8.705054	161.595047	2228.305664	175.472198	234.087051	983.675720	-0.342254	-0.294564	1.821556	1.884842	-0.513369	0.651058	0.537046	-2.615676	-2.487624	4.619804	4.233167	5.308126	-5.576128	-2.581536	-4.365026	-6.566184	-1.339482	-2.204539	-1.126629	-2.361625	-2.628877	-0.362707	0.566276	-0.091297	0.440614	1.080030	-0.755546	-0.407814	-1.183697	-0.848132	-0.313204	-2.134031	8.850863	2.711628	4.811763	8.958502	-0.366285	3.318763	0.329840	1.538568	4.348969	-1.816757	6.563427	1.253030	2.842146	7.301911	-0.718124	2.376652	-0.260213	0.996117	2.531418	-1.418708	2.053362	1.714995	5.849802	5.440816	1.884865	1.741266	2.375170	0.947722	3.170853	1.950660	286937	0
834854	168	0.596468	0.002326	0.002032	0.005509	0.007654	0.003616	-1	-1.436731	-2.133809	4.648743	2.969781	0.907301	0.450538	1.191757	-0.143514	9.682172	8.077853	18.023613	10.920301	25.319948	17.384563	-1.470224	-1.229787	-3.497435	-1.578663	-4.015155	-3.872986	-5.276647	-4.489394	4.190367	3.187598	10.996456	9.465394	19.159149	12.718694	31.042587	19.925367	-0.962800	-0.650292	4.028402	4.142321	0.611631	1.429632	2.867671	-2.501698	-0.290131	-0.013775	-0.325107	-0.721573	-0.889294	-0.459434	-0.440759	-0.153293	7.292922	6.501546	23.773369	1444.840088	87.525017	303.269165	1372.886841	1.325848	1.173081	-0.216417	-0.236258	1.529365	2.043775	1.743800	3.031896	2.989350	-0.761621	-2.308117	-1.328938	0.220087	0.031546	0.051950	0.083291	0.902025	0.278071	0.052404	0.093597	0.077576	1.635452	-0.001564	-1.001347	-1.168230	-0.163837	-0.383313	0.311700	-1.112664	-0.479528	0.244936	0.945099	-0.420525	-1.712330	-0.716374	-0.219099	-0.496815	-0.363410	-0.959098	-0.536406	0.115152	0.687828	-1.948212	-1.659458	-1.141388	-1.382763	-0.595652	-1.113544	-1.988687	-1.147521	-0.065230	0.619337	0.335134	0.268782	7.619550	7.574471	0.335153	0.244880	0.339170	0.232385	8.562437	7.662621	834854	1

データセットには極端な外れ値があると推測できます。また、一部の外れ値には、特徴間の高い多重共線性の結果である隣接列からの大きな値が付随していることにも言及する価値があります。

次に、特徴量データの99.9％を超えるこれらの極端な外れ値を取り除いてみましょう。

データセットのループ中にデータが削除されないようにするために、すべての特徴量に対して99.9％をマークするためのリストを作成します。また、後で調べるために「.1％のマークを使用して」という負の外れ値のリストを作成します。

n999 = [ np.percentile(df[i],99.9) for i in featstr[1:]]
n001 = [ np.percentile(df[i],.1) for i in featstr[1:]]

for i, j in enumerate(featstr[1:]):
    df = df[df[j] < n999[i]]
    gc.collect()

どの程度のデータを削除したか見ましょう。

str(round(((org_len - len(df))/org_len)*100,2))+'%'

'4.01%'

fig = px.bar(df[featstr].max(), title='Features Max Values')
fig.layout.xaxis.tickangle = 300
fig.update_traces(showlegend = False)
fig.layout.xaxis. dtick = 5
fig.layout.xaxis.title = ''
fig.layout.yaxis.title = ''
fig.show()

次に、カスタマイズされた.1％と99.9％のひげを使用して箱ひげ図グリッドを再度作成します。

fig = plt.figure(figsize=(20,80))
fig.suptitle('Features Box plot with 0.1% 99.9% whiskers',fontsize=22, y=.89)
grid =  gridspec.GridSpec(29,4,figure=fig,hspace=.5,wspace=.05)
counter = 0
for i in range(29):
    for j in range(4):
        subf = fig.add_subplot(grid[i, j]);
        sns.boxplot(x= df[featstr[counter]],saturation=.5,color= 'blue', ax= subf,width=.5,whis=(.1,99.9));
        subf.set_xlabel('')
        subf.set_title('{}'.format(featstr[counter]),fontsize=16)
        counter += 1
        gc.collect();

正の外れ値のみを削除したため、データセットにはまだ多くの明らかな外れ値、特に負の値があります。そのため、以前は対称的な外れ値を持っていた特徴量に、ある種の左に歪んだ外れ値があります。
両方の箱ひげ図を比較すると（正の外れ値を削除する前と後）、以前は対称的な外れ値を持っていたfeature_3からfeature_40までの特徴量が、トリミング後に極端な負の外れ値を持つように変更されていることがわかります。

手動でこれらの特徴量をトリミングします。

for i,j in zip(featstr[1:][2:34],n001[2:34]):
    df = df[df[i] > j]
    gc.collect();

次のコードは、これまでのクリーニングのために失われた観測値の比率を示しています。

str(round(((org_len - len(df))/org_len)*100,2))+'%'

'4.89%'

外れ値を削除し、resp値を考慮した後の特徴量のkde(カーネル密度推定)とヒストグラムを見ましょう。

fig = plt.figure(figsize=(20,80))
fig.suptitle('KDE plot of Features',fontsize=24,transform =fig.transFigure, y=.89)
grid =  gridspec.GridSpec(29,4,figure=fig,hspace=.5,wspace=.01)
counter = 0
for i in range(29):
    for j in range(4):
        subf = fig.add_subplot(grid[i, j]);
        sns.distplot(df[df.action==0][featstr[counter]],bins= 100,label='Negative',
                     color='darkorange', kde_kws={'linewidth':4},ax=subf)
        sns.distplot(df[df.action!=0][featstr[counter]],bins= 100,label='Positive',
                     color='blue', kde_kws={'alpha':.9,'linewidth':2},hist_kws={'alpha':.3},ax=subf)
        subf.axvline(np.percentile(df[featstr[counter]],99.5),color= 'darkblue', label='99.5%', linestyle=':',linewidth=2)
        subf.axvline(np.percentile(df[featstr[counter]],.5),color= 'red', label='0.5%', linestyle=':',linewidth=2)
        subf.legend().set_visible(False)
        subf.set_xlabel('')
        subf.set_title('{}'.format(featstr[counter]),fontsize=16)
        kurt=round(df[featstr[counter]].kurt(),2)
        skew=round(df[featstr[counter]].skew(),2)
        subf.text(.6,.92,'Kurt = {:.2f}\nSkew = {:.2f}'.format(kurt ,skew),
         transform=subf.transAxes, verticalalignment='top',bbox=props,fontsize=10)
        counter += 1
        gc.collect();
handles, labels = subf.get_legend_handles_labels()
fig.legend(handles, labels,ncol=4, bbox_to_anchor=(0.86, 0.893),fontsize=10,
           title= 'Resp',title_fontsize=14,bbox_transform =fig.transFigure);

特徴量の分布にresp値を追加すると、次のようになります。
- 特徴のヒストグラムは、外れ値がはるかに少なくなり、より正式な分布になりました。
- また、1、2、85、87、88、91などの一部の機能には多くの負の外れ値があることがわかります。
- 49、50、51、55、56、57、58、59などの一部の機能には、まだ多くの正の外れ値があります。
- 特徴量の分布は、resp値の影響を受けません。

特徴量とrespの相関

まず、respと各特徴量の関係についてpandasの相関系列を作成します。

respcorr =  pd.Series([ df.resp.corr(df[i]) for i in featstr],index=featstr)

fig = px.bar(respcorr,color = respcorr, color_continuous_scale=['red','blue'], title= 'Features Correlation with Resp')
fig.layout.xaxis.tickangle = 300
fig.layout.xaxis. dtick = 5
fig.layout.xaxis.title = ''
fig.layout.yaxis.title = 'pearson correlation'
fig.update(layout_coloraxis_showscale=False)
fig.show();

特徴量がrespとそこまで高い相関がないことがわかります。

特徴量と重み(weight)

ここでも、pandasの一連の重みと特徴量の相関を見ますが、0より大きい重みのみを考慮します。

wecorr = pd.Series([df[df.weight != 0].weight.corr(df[df.weight != 0][i]) for i in featstr],index=featstr)

wecorr.head(10)

feature_0    -0.004811
feature_1    -0.039106
feature_2     0.061847
feature_3     0.000140
feature_4     0.005604
feature_5    -0.001310
feature_6     0.001633
feature_9    -0.096143
feature_10    0.012905
feature_11   -0.041127
dtype: float64

fig = px.bar(wecorr,title= 'Features Correlation with Weight (not including zero weights)')
fig.layout.xaxis.tickangle = 300
fig.layout.xaxis. dtick = 5
fig.layout.xaxis.title = ''
fig.layout.yaxis.title = 'pearson correlation'
fig.update(layout_coloraxis_showscale=False)
fig.update_layout(showlegend=False)
fig.show()

散布図を使用して、相関係数の最大値および最小値を調べます。

まずは、重みに対して最も高い相関を持つ特徴量(feature_51)について見ます。

fig = plt.figure(figsize=(8,6))
sns.scatterplot(df[df.weight != 0].weight,df[df.weight != 0].feature_51, color = 'darkblue', alpha=.3)
plt.xlabel('Weight',fontsize=14)
plt.ylabel('Featre_51',fontsize=14)
plt.title('Feature_51 vs Weight\nCorrelation = {}%'.format(round(df[df.weight != 0].weight.corr(df[df.weight != 0].feature_51),4)*100),fontsize=16);

重み(weight)はfeature_51と散布図からも高い相関があることがわかります。

重みに最も低い相関を持つ特徴量(feature_126)について見ましょう。

fig = plt.figure(figsize=(8,6))
sns.scatterplot(df[df.weight != 0].weight,df[df.weight != 0].feature_126, color = 'darkblue', alpha=.3)
plt.xlabel('Weight',fontsize=14)
plt.ylabel('Featre_126',fontsize=14)
plt.title('Feature_126 vs Weight\nCorrelation{}%'.format(round(df[df.weight != 0].weight.corr(df[df.weight != 0].feature_126),4)*100),fontsize=16);

重みとfeature_126の間にはある種の負の相関関係がありますが、関係は弱いようです

Feature 0

feature_0の一意の値を見つけましょう。

plt.figure(figsize=(7,5)) 
df.feature_0.value_counts().plot.bar(color='darkblue',alpha=.6,width=.5)
plt.title('Feature_0',fontsize=18) 
plt.xticks(rotation=0,fontsize=14);

Feature_0はバイナリな特徴量のようです

Respを考慮する。

plt.figure(figsize=(8,6)) 
sns.countplot(data=df, x='feature_0', hue='action',palette='viridis')
plt.xticks(fontsize=14)
plt.yticks(fontsize=13)
plt.xlabel('Feature 0',fontsize=12)
plt.title('Feature 0 and Resp', fontsize=18)
plt.ylabel('')
plt.xlim(-1,2)
h, l = plt.gca().get_legend_handles_labels()
plt.legend(h,['Negative','Positive'],ncol=1, fontsize=12, loc=3,title= 'Resp',title_fontsize=14);

feature_0とrespが負または正である間に明確な関係がないようです