[Ubiquant Market Prediction]特徴量エンジニアリングのアイデア

こんにちは、tmrtjと申します。

私は現在株式会社音圧爆上げくんという会社に所属しており、プロKagglerとして活動しています。


私は現在プロKagglerとして、業務の一環で、現在Kaggleで行われているコンペの一つである[Ubiquant Market Prediction](https://www.kaggle.com/c/ubiquant-market-prediction)に参加しています。



また、Kaggle上にCodeとして同じ内容の記事を英語で投稿しています。
リンクは[こちら](https://www.kaggle.com/tmrtj9999/feature-engineering-idea)です。

#コンペについて
Ubiquant Market Prediction は匿名化されたデータセットが与えられており、それらを用いてtargetを予測するコンペです。
<br>
特徴量としては300個の匿名化されたデータが与えられていますが、生の資産価格データが与えられているわけではないので、通常の金融時系列データ分析で行うような特徴量エンジニアリングが困難になっています。
<br>
そこで、今回は、新たな特徴量を作るアイディアを一つ紹介します。

#考え方
EDAとして多くの人が特徴量の分布を確認することを行っていると思います。例えば、このようにヒストグラム化できます。これは'f14'のヒストグラムです。


![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2465043/87afe00b-b526-55cf-ba04-0ea7f8780463.png)


これらの特徴量の分布を眺めていると、分布がきれいに左右対称になっているものもあれば、非対称になっているものもあります。実際の計算式は分かりませんが、特徴量の計算の仕方でこのような特徴が現れるのだと思われます。そこで、これらの特徴量の数値が、どれくらい真ん中から外れているかを特徴量に組み込むことを考えました。具体的な計算式としては、

（特徴量の値　- 特徴量の中央値)の絶対値

として計算します。
<br>
以下で実装していきます。


#実装

まずは下準備をします。
trainデータを読み込んでDataFrameにします。
```reduce_mem_usage(df)```　はメモリを節約できるような変換をかける関数です。


```
import os
import pandas as pd
import numpy as np
from sklearn.model_selection import train_test_split
import matplotlib.pyplot as plt
%matplotlib inline

import warnings
warnings.filterwarnings('ignore')

import lightgbm as lgbm
from lightgbm import *


df = pd.read_parquet('../input/ubiquant-parquet/train_low_mem.parquet')

def reduce_mem_usage(df):
  
    start_mem = df.memory_usage().sum() / 1024**2
    print('Memory usage of dataframe is {:.2f} MB'.format(start_mem))
    
    for col in df.columns:
        col_type = df[col].dtype
        
        if col_type != object:
            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 < np.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)
        else:
            df[col] = df[col].astype('category')

    end_mem = df.memory_usage().sum() / 1024**2
    
    print('Memory usage after optimization is: {:.2f} MB'.format(end_mem))
    print('Decreased by {:.1f}%'.format(100 * (start_mem - end_mem) / start_mem))
     
    return df

df = reduce_mem_usage(df)

```


<br>
'f_146'はヒストグラム化すると図のようになっています。

```
df['f_146'].hist(bins = 100)
```


![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2465043/32e2fdb7-5129-e9a1-fdb0-e015cf0da3a9.png)

<br>
<br>
これは見ての通り綺麗な左右対称の分布になっています。
次に、特徴量を加工します。
今回は'f146'を加工したものが効いているのかを確認しますが、せっかくなので一斉に中央値からの中央値からの差を計算します。



```
features = [f'f_{i}' for i in range(300)]
df_median = df[features].head(30000).median()
for i in range(300):
    df[f'f_median_{i}'] = abs(df[f'f_{i}']-df_median[f'f_{i}'])
```


<br>
<br>

実際にこの特徴量が効いているのかを確認するために、LightGBMを使用し、特徴量重要度を確認します。


中央値計算に使用した行は学習用データに含まれないようにしながらtrain, val, testデータに分割します。

```
from sklearn.model_selection import KFold, train_test_split


df = df.tail(3101410)
features = [f'f_{i}' for i in range(300)] + [f'f_median_146']
target = 'target'

df_features = df[features]

X_train, X, Y_train, Y = train_test_split(df_features, df[target], train_size=0.6, shuffle=False)

df = [[]]
df_features = [[]]

X_val, X_test, Y_val, Y_test = train_test_split(X, Y, train_size=0.5, shuffle=False)
```

<br>
LightGBMにかけます。
<br>

```
import warnings
import numpy as np
import lightgbm as lgb
from scipy.stats import pearsonr

warnings.simplefilter('ignore')

lgb_train = lgb.Dataset(X_train, Y_train)
lgb_eval = lgb.Dataset(X_val, Y_val, reference=lgb_train)

params = {#'verbosity' : -1,
            #'extra_trees': True,
            #'boosting': "rf",
            'seed': 1,
            'objective': "regression",
            'learning_rate': 0.02,
            'max_depth' : 6,
            'num_leaves': 64,
            'feature_fraction': 0.1, #0.05の方が精度良かった
            'bagging_fraction': 0.1,
            #'bagging_freq': 2
          }
        
        
gbm = lgb.train(params,
                lgb_train,
                num_boost_round=19450815,
                valid_sets=lgb_eval,
                #verbose_eval=False,
                early_stopping_rounds=10,
                )


Y_pred = gbm.predict(X_test, num_iteration=gbm.best_iteration)

score_tuple = pearsonr(Y_test, Y_pred)
score = score_tuple[0]
print(f"Validation Pearsonr score : {score:.4f}")
```

**結果**

```
Validation Pearsonr score : 0.1193
```

<br>
<br>
特徴量重要度を表示します。


```
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.ensemble import RandomForestClassifier
from sklearn.datasets import load_iris
from sklearn.model_selection import train_test_split
from sklearn.metrics import classification_report

#特徴量の重要度
feature = gbm.feature_importance(importance_type='gain')

#特徴量の重要度を上から順に出力する
f = pd.DataFrame({'number': range(0, len(feature)),
             'feature': feature[:]})
f2 = f.sort_values('feature',ascending=False)

#特徴量の名前
label = X_train.columns[0:]

#特徴量の重要度順（降順）
indices = np.argsort(feature)[::-1]

for i in range(len(feature)):
    print(str(i + 1) + "   " + str(label[indices[i]]) + "   " + str(feature[indices[i]]))

```
<br>
<br>
**結果**

```
// ... 略
266   f_262   547.4649114608765
267   f_median_146   542.8843946456909
268   f_24   531.7534027099609
269   f_124   530.3465976715088
// ... 略

```







特徴量重要度を確認すると、確かに"f_median_146"が効いていることがわかります。

他にも、分布が左右対称の特徴量はあるので同様のやり方で使える特徴量の加工できるものはあると思われます。


#さいごに
株式会社音圧爆上げくんではプロKagglerを募集しています。
ご興味のある方はぜひ以下のリンクをご覧ください。
[Wantedlyリンク](https://www.wantedly.com/projects/748285)