【kaggleで機械学習勉強・第一回】Logistic Regression【kaggle,python,ロジスティック回帰】

今回はchris deotteさんのカーネル
Logistic Regression - 0.800　を紹介しつつ、本文とコードを翻訳していく。
chris deotteさんはコンペに対して複数のアプローチをしていて、
私自身も多くの手法を知っていきたいと思っていたのでその全てを追いかけていく。(全5回くらい)
そのため今回は第一回と題している。

コンペティションのテーマ

今回のコンペのinstant gratificationは解釈が難しい部分が多く正確に翻訳できた自信がない。
翻訳に失敗している場合は指摘していただきたい。

instant gratificationで使用するデータは、
流れ着いた瓶の中にUSBが入っていて変な名前の列名と何の値かわからない値が入っていた。
というような内容。

さらに、通常のコンペではtestのデータを予測するものだが、
今回はtestのデータも欠けており、kaggle側のみが正解と欠けた残りのデータを持っている。
欠けた部分を補完したデータはプライベートデータセットと呼ばれる。
この欠けた部分も考慮しながらモデルを作成しなくてはならない。
予測したいデータの構造が不明な状態でどう予測していくのか？

早速

ロジスティック回帰を使った結果、リーダーボードの結果が0.800だった。

結論から言って、変数「wheezy-copper-turtle-magic」というものが他の変数と関係を持っている。

単純なモデルながら、良く機能している。

変数間を独立とした場合のロジスティックス回帰ではスコア0.5であったが、
変数間の関係を考慮した場合のロジスティックス回帰では0.800のスコアでした。

import numpy as np, pandas as pd, os
from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import StratifiedKFold
from sklearn.metrics import roc_auc_score

train = pd.read_csv('C:/train.csv')
test = pd.read_csv('C:/test.csv')

train.head()

id                       muggy-smalt-axolotl-pembus     dorky-peach-sheepdog-ordinal・・・
707b395ecdcbb4dc2eabea00e4d1b179    -2.070654   1.018160・・・
5880c03c6582a7b42248668e56b4bdec    -0.491702   0.082645・・・

変数の名前はコンペの説明にあった通り不明な名前だった。
258列のうちIDとtargetを除く256の変数がある。
ターゲットとの相関を見てみると、0.04を超える相関を持った変数は無い

なので、変数間を独立のものとしてロジスティックス回帰をターゲットに対して行うと、スコアは約0.5だった。

wheezy-copper-turtle-magicが特別な変数として扱われるけれど、これはこの変数だけが質的変数だから。
(このあたりは本カーネルでは言葉の上でしか表現されていない。EDAで見つけるのは他の人のカーネルを参照すべし。)

全変数をplotしてみてもわかるが、EDAのついでにpandas-profilingで特定するのが楽かと思ってやってみた。

import pandas as pd
import pandas_profiling as pdp

h20 = train.head(20)
pdp.ProfileReport(h20)

データが重いので頭部分から20行しかやってないけれど、平均値以外の統計量が整数で、質的変数でありそうな感じがする。

cols = [c for c in train.columns if c not in ['id', 'target']]
oof = np.zeros(len(train))
skf = StratifiedKFold(n_splits=5, random_state=42)

for train_index, test_index in skf.split(train.iloc[:,1:-1], train['target']):
 clf = LogisticRegression(solver='liblinear',penalty='l2',C=1.0)
 clf.fit(train.loc[train_index][cols],train.loc[train_index]['target'])
 oof[test_index] = clf.predict_proba(train.loc[test_index][cols])[:,1]

auc = roc_auc_score(train['target'],oof)
print('LR without interactions scores CV =',round(auc,5))

LR without interactions scores CV = 0.52994

コードがきになったので。

colsはカラムの名前
oofという0で埋められたハコをつくり、
skfに交差検証用の層状k分割を指示する
(データを訓練とテスト用に分割している)

for文で、
・train.iloc[:,1:-1]　idとターゲットを除いたもの
・train['target']　ターゲット
これをskf.splitで分割して、
train_index, test_index　にそれぞれ入れて

sklearnのLRでlibkinearで最適化を選択し、
trainに適応させて
offに予測値を入力して、
予測値と実際の値のaucを比較して、その値を精度としている。

参考サイト

https://data-science.gr.jp/implementation/iml_sklearn_logistic_regression.html
https://medium.com/music-and-technology/scikit-learn%E3%81%AEstratifiedkfold%E3%81%A7%E4%BA%A4%E5%B7%AE%E6%A4%9C%E5%AE%9A-fcb3c9f4572b

LogisticRegressionのsolverが何をやっているかわからないんです。って質問
解説3の項目を参照
https://stackoverflow.com/questions/38640109/logistic-regression-python-solvers-defintions

関係を考慮して512の別々モデルを作成

wheezy-copper-turtle-magicという変数が質的変数として存在しており、
この変数に対応して他の変数は分布が異なるようになっているのでは？
ということを考えている。
wheezy-copper-turtle-magicは512個の整数で成り立っているので、
wheezy-copper-turtle-magic==1の時のモデル
wheezy-copper-turtle-magic==2の時のモデル
・・・
といった具合に作る。

結論から言うとこのモデルの作成方法ではスコアが0.8まで上がる。

cols.remove('wheezy-copper-turtle-magic')
interactions = np.zeros((512,255))
oof = np.zeros(len(train))
preds = np.zeros(len(test))

# BUILD 512 SEPARATE MODELS
for i in range(512):
    # ONLY TRAIN WITH DATA WHERE WHEEZY EQUALS I
    train2 = train[train['wheezy-copper-turtle-magic']==i]
    test2 = test[test['wheezy-copper-turtle-magic']==i]
    idx1 = train2.index; idx2 = test2.index
    train2.reset_index(drop=True,inplace=True)
    test2.reset_index(drop=True,inplace=True)

    skf = StratifiedKFold(n_splits=25, random_state=42)
    for train_index, test_index in skf.split(train2.iloc[:,1:-1], train2['target']):
        # LOGISTIC REGRESSION MODEL
        clf = LogisticRegression(solver='liblinear',penalty='l1',C=0.05)
        clf.fit(train2.loc[train_index][cols],train2.loc[train_index]['target'])
        oof[idx1[test_index]] = clf.predict_proba(train2.loc[test_index][cols])[:,1]
        preds[idx2] += clf.predict_proba(test2[cols])[:,1] / 25.0
        # RECORD INTERACTIONS
        for j in range(255):
            if clf.coef_[0][j]>0: interactions[i,j] = 1
            elif clf.coef_[0][j]<0: interactions[i,j] = -1
    #if i%25==0: print(i)

# PRINT CV AUC
auc = roc_auc_score(train['target'],oof)
print('LR with interactions scores CV =',round(auc,5))

ここではそんなに複雑なことはやっていなくて、
データをwheezy-copper-turtle-magicの値ごとに切り分けて取り出して、
その時の各変数をskfで切り出して、ロジスティック回帰にかけているだけ。

LR with interactions scores CV = 0.80549

512個のモデルをつくったところスコアは0.80まで上がった。

sub = pd.read_csv('sample_submission.csv')
sub['target'] = preds
sub.to_csv('submission.csv',index=False)

submitできるようにして終了。

まとめ

wheezy-copper-turtle-magicが、targetを予測する際に他の変数との関係をもっている変数であることが分かった。
それと(ロジスティック回帰のような)シンプル(単純)なモデルでも高いスコアを出すことが出来ると分かった。

今回のカーネルの考え方は排他的論理和(XOR?ビット演算?)の問題に似ている。
二つの変数があって、一つのターゲットがある。(x1,x2,y)とする。
x1もx2も同じビットなら0を返し、異なれば1を返すので
(0,0,0), (1,0,1), (0,1,1), (1,1,0)
x1やx2はyに相関関係を持たないことは注意しておいてください。
x1とx2の間にも相関を持っていません。
しかしx1とx2が同じでなければyは0にならないし、異なっていなければ1にならないという関係性はもっています。
(こんな翻訳であってるのか？)

関係性

wheezy-copper-turtle-magicとほかの変数の関係について記述します。
各変数をまとめてモデルを作るだけではtargetを予測できませんが、
wheezy-copper-turtle-magicが「特定の値」であるとき、ほかの変数を使ってtargetを予測できるんです。
たとえば、wheezy-copper-turtle-magic = 0の場合、
zippy-harlequin-otter-grandmasterはtargetと正の相関があります。
wheezy-copper-turtle-magic = 0の場合、hasty-puce-fowl-fepidはターゲットと負の相関があります。

import seaborn as sns
import matplotlib.pyplot as plt
import warnings
warnings.filterwarnings("ignore")

plt.figure(figsize=(15,5))

# PLOT ALL ZIPPY
plt.subplot(1,2,1)
sns.distplot(train[ (train['target']==0) ]['zippy-harlequin-otter-grandmaster'], label = 't=0')
sns.distplot(train[ (train['target']==1) ]['zippy-harlequin-otter-grandmaster'], label = 't=1')
plt.title("Without interaction, zippy has no correlation \n (showing all rows)")
plt.xlim((-5,5))
plt.legend()

# PLOT ZIPPY WHERE WHEEZY-MAGIC=0
plt.subplot(1,2,2)
sns.distplot(train[ (train['wheezy-copper-turtle-magic']==0) & (train['target']==0) ]
             ['zippy-harlequin-otter-grandmaster'], label = 't=0')
sns.distplot(train[ (train['wheezy-copper-turtle-magic']==0) & (train['target']==1) ]
             ['zippy-harlequin-otter-grandmaster'], label = 't=1')
plt.title("With interaction, zippy has postive correlation \n (only showing rows where wheezy-copper-turtle-magic=0)")
plt.legend()

plt.show()

やっぱseaborn使うと奇麗になるな。
上の図が正の相関を持っている図。

plt.figure(figsize=(15,5))

# PLOT ALL HASTY
plt.subplot(1,2,1)
sns.distplot(train[ (train['target']==0) ]['hasty-puce-fowl-fepid'], label = 't=0')
sns.distplot(train[ (train['target']==1) ]['hasty-puce-fowl-fepid'], label = 't=1')
plt.title("Without interaction, hasty has no correlation \n (showing all rows)")
plt.xlim((-5,5))
plt.legend()

# PLOT HASTY WHERE WHEEZY-MAGIC=0
plt.subplot(1,2,2)
sns.distplot(train[ (train['wheezy-copper-turtle-magic']==0) & (train['target']==0) ]
             ['hasty-puce-fowl-fepid'], label = 't=0')
sns.distplot(train[ (train['wheezy-copper-turtle-magic']==0) & (train['target']==1) ]
             ['hasty-puce-fowl-fepid'], label = 't=1')
plt.title("With interaction, hasty has negative correlation \n (only showing rows where wheezy-copper-turtle-magic=0)")
plt.legend()

plt.show()

負の相関

そんな感じで他の変数にも関係性があって、512のモデルを作った時、256個の変数は以下の図のような相関関係で示される。
青が負の相関で黄色が正の相関である。

# PLOT INTERACTIONS WITH WHEEZY-MAGIC
plt.figure(figsize=(15,8))
plt.matshow(interactions.transpose(),fignum=1)
plt.title("Variable Interactions with wheezy-copper-turtle-magic \n \
    Yellow = combines to create positive correlation with target. Blue = negative correlation",fontsize=16)
plt.xlabel('values of wheezy-copper-turtle-magic')
plt.ylabel('columns of train')
plt.show()

x1変数が1,x2変数が0,x3変数が1......というようなパターン性からtargetが1か0かをXORのような仕組みで予測した。
ってことを言いたかったのかな？
その1,0の数値というのが相関関係であって、
wheezy-copper-turtle-magic変数が0~511までの値をとるとき、
その他のすべての変数の相関関係が異なってくるから高い精度で判別できた。
ということだと私は翻訳していて思いました。

ここからはコメントコーナーを翻訳

Q:なんで最初のモデルにL2を使って、二番目のモデルではL1をつかってるの？
(ロジスティック回帰のオプションで罰則を選択する部分の話と思われる)

A:じつは最初のモデルでも二番目のモデルでもL1,L2どちらも試していて、結果の良いペナルティを選択した。比較できるように本文中に結果を書いておくべきでした。

Q:今回はtargetと相関を持たない変数達だったみたいですが、もし一つもしくは複数がtargetと相関をもっていたらどうしますか?

A:今回のカーネルでやった方法はもちろん試しますし、相関を考えたモデルも作って、モデル同士をくっつけて予測モデルを構築すると思います。

Q:結果の判断は正しいの？

A:ROCやAUCは判断基準としていいとおもう。
AUCの尺度は　正しくペアリングされたもの/全体のペア数　を表しているから。
確率を計算するには、予測を25で割る必要があります。

preds[idx2] += clf.predict_proba(test2[cols])[:,1] / 25.0

でも、vladislavさんが指摘しているように、判断基準はAUCであり、AUCはすべての予測値の順序関係を考慮します。
混乱を防ぐためにもカーネルの中でも25で割ったものを示すべきでした。

Q:なぜモデルの数でOOF配列のすべての要素を分割しなかったのですか？

A:OOFのなかで更新させているのです。
要素500個でモデル1をつくって、
また別の要素500個で書き換えて、それをモデル2として。のように。

Q:関係性を説明するときの分布曲線から正負の相関はどうみたら？

A:赤色のヒストグラムがx軸の正の方向にシフトしていってることから正の相関
　hastyの変数は左にシフトしてるので負の相関
　zippyの変数はターゲットが1の時正の相関で、0のとき負の相関。

Q:magic変数(例の質的変数のこと)でモデルを区切ろうってどうして思いついたの？

A:たとえばtitanicでも性別でセグメント分けするでしょ。
　性別が予測にどう影響するかはEDAの段階で考えなきゃいけない部分でしょう。

Q:skfを25にした理由は?

A:magik変数=kの時、トレーニングとして食わせる値は500個になります。
foldingを5にしたら一度のトレーニングは400個をつかって100個は使われなくなります。
foldingを25にすると480個を使って学習してくれます。
256個も変数があるので、トレーニングに使用するものがわずかに異なるだけでモデルや結果が変わってきます。
試しにfoldingを25以下にしたらスコアがさがりました。
folding=15; cv=0.804
folding=11; cv=0.803
folding=7; cv=0.802

Q:今回ロジスティック回帰のハイパーパラメーター(事前に与えられるパラメーター)のグリッドサーチとかどんなことをやりましたか？
A:グリッドサーチはしてない。512のモデルに同じパラメーターを使いました。グリッドサーチしたらさらに良いモデルができるかもしれませんね。

以上

つづく