えんじにあのじゆうちょう

勉強したことを中心にアウトプットしていきます。

SIGNATEの【第2回_Beginner限定コンペ】健診データによる肝疾患判定をやってみた

はじめに

SIGNATEさんは去年こちらの講座でお世話になった以来放置していたので、相変わらずBeginnerだったのですが、Beginner限定コンペであるラインを超えればIntermediateに昇格ということだったのでやってみることにしました。
せっかくなのでそのログとしてこの記事を残します。

取り組んだこと

実際にやったこと

いちいち順は追いませんが、データを見て特徴量作って、モデルを当てはめてチューニングしたっていう感じです。
最終的には9/4昼頃の段階で0.9271133というスコアを出し無事クリアしたという感じです。

ちなみにモデルは最初XGBoostで頑張っていたのですが、0.918くらいで止まってしまったので全く同じデータのままlightgbmに切り替えたらサクッと超えてしまったという感じです。

特徴量生成

最初は真面目に可視化していました。例えばこれ。

headers = train_df.drop(columns=["Gender", 'disease']).columns
fig, ax = plt.subplots(headers.shape[0], headers.shape[0], figsize=(20, 20))
for i, h1 in enumerate(headers):
    for j, h2 in enumerate(headers):
        ax[i, j].scatter(disease_df[h1], disease_df[h2], color="red")
        ax[i, j].scatter(not_disease_df[h1], not_disease_df[h2], color="blue")
plt.show()

f:id:marufeuillex:20200905134316p:plain
(途中まで)

赤丸がdisease=1の点です。
うーん、だいたい値が大きくなると赤丸増えるなーと思って気づきました。
そうだ。これは検診データだ。これらが正確に何かはわからないけど、検診データということは基本的に正常範囲があるはず、と。
そんなわけで、正常範囲を調べたりしながら、特徴量を作ることにしました。

また、次元削減した特徴量も効くことがある、ということで、いくつか次元削減を試しました。

まずはt-SNE。rand_seedがいい感じのやつを探そう!と思い100くらいまで表示してみましたが、結局何も考えずにゼロを選んでいたりします。

header = ["latent{}".format(i) for i in range(2)]
sc = StandardScaler()
sc.fit(train.drop(columns=["disease"]))
fig, ax = plt.subplots(10, 10, figsize=(20, 20))
disease_idx = train["disease"] == 1
not_disease_idx = train["disease"] != 1
for i in range(100):
    latent_features = pd.DataFrame(
        bhtsne.tsne(
            sc.transform(
                train.drop(columns=["disease"])
            ),
            dimensions=2,
            rand_seed=i
        ),
        columns=header
    )
    disease_latent = latent_features[disease_idx]
    not_disease_latent = latent_features[not_disease_idx]
    ax[i//10, i%10].scatter(disease_latent["latent0"], disease_latent["latent1"], color="red")
    ax[i//10, i%10].scatter(not_disease_latent["latent0"], not_disease_latent["latent1"], color="blue")
plt.show()

f:id:marufeuillex:20200905134727p:plain

ちなみにUMAPも試してみたのですが、残念ながらこれよりパッと見悪かったので不採用としています。

あと、PCAも試してみました。

pca = PCA(n_components=2)
sc = StandardScaler()
feature_sc = sc.fit_transform(train_bst_df.drop(columns=["disease"]))
latent = pca.fit_transform(feature_sc)
disease_idx = train_bst_df["disease"] == 1
not_disease_idx = train_bst_df["disease"] != 1
plt.scatter(latent[disease_idx][:, 0], latent[disease_idx][:, 1], color="red")
plt.scatter(latent[not_disease_idx][:, 0], latent[not_disease_idx][:, 1], color="blue")

f:id:marufeuillex:20200905134844p:plain
微妙な気もしますが、木で分割できそうな気もするので、一応これも採用しています。

そんなこんなで、結局いかのように特徴量をつけています。

def add_features(train, test):
    train_df = train.copy()
    test_df = test.copy()
    # add feature
    disease = np.zeros(9)
    total = np.zeros(9)

    # 年齢ごとにdiseaseの率に偏りがあったため、それを表す特徴を作成
    for age, dis in train[["Age", "disease"]].to_numpy():
        b = age // 10
        disease[b] += dis
        total[b] += 1
    # 0除算を避ける目的
    total[total == 0] = 1
    r = np.nan_to_num(disease / total)
    
    df = pd.concat([train_df, test_df]).reset_index(drop=True)
    
    df["Age_Feature"] = df["Age"].map(lambda age: r[age//10])
    
    
    # Genderは文字列なので、0,1に変換
    ohe = OneHotEncoder(categories="auto", sparse=False, dtype=np.int, drop='first')
    df["Gender_OneHot"] = ohe.fit_transform(df[["Gender"]])
    df = df.drop(columns=["Gender"])
    
    # 一般的な指標内かどうか
    # 指標の上限 or 下限からどの程度離れているか
    df["AG_ratio_error"] = 0
    df["AG_ratio_outlier"] = 0
    df.loc[(df["AG_ratio"] < 1.1) | (2.0 < df["AG_ratio"]), "AG_ratio_error"] = 1
    df.loc[(df["AG_ratio"] < 1.1), "AG_ratio_outlier"] = df["AG_ratio"] - 1.1
    df.loc[(2.0 < df["AG_ratio"]), "AG_ratio_outlier"] = df["AG_ratio"] - 2.0
    
    df["T_Bil_error"] = 0
    df["T_Bil_outlier"] = 0
    df.loc[(df["T_Bil"] < 0.2) | (1.2 < df["T_Bil"]), "T_Bil_error"] = 1
    df.loc[(df["T_Bil"] < 0.2), "T_Bil_outlier"] = df["T_Bil"] - 0.2
    df.loc[(1.2 < df["T_Bil"]), "T_Bil_outlier"] = df["T_Bil"] - 1.2
    
    df["ALP_error"] = 0
    df["ALP_outlier"] = 0
    df.loc[(df["ALP"] < 50) | (350 < df["ALP"]), "ALP_error"] = 1
    df.loc[(df["ALP"] < 50), "ALP_outlier"] = df["ALP"] - 50
    df.loc[(350 < df["ALP"]), "ALP_outlier"] = df["ALP"] - 350
    
    df["ALT_GPT_error"] = 0
    df["ALT_GPT_outlier"] = 0
    df.loc[(df["ALT_GPT"] < 5) | (40 < df["ALT_GPT"]), "ALT_GPT_error"] = 1
    df.loc[(df["ALT_GPT"] < 5), "ALT_GPT_outlier"] = df["ALT_GPT"] - 5
    df.loc[(40 < df["ALT_GPT"]), "ALT_GPT_outlier"] = df["ALT_GPT"] - 40
    
    df["AST_GOT_error"] = 0
    df["AST_GOT_outlier"] = 0
    df.loc[(df["AST_GOT"] < 10) | (40 < df["AST_GOT"]), "AST_GOT_error"] = 1
    df.loc[(df["AST_GOT"] < 10), "AST_GOT_outlier"] = df["AST_GOT"] - 10
    df.loc[(40 < df["AST_GOT"]), "AST_GOT_outlier"] = df["AST_GOT"] - 40
    
    df["TP_error"] = 0
    df["TP_outlier"] = 0
    df.loc[(df["TP"] < 6.7) | (8.3 < df["TP"]), "TP_error"] = 1
    df.loc[(df["TP"] < 6.7), "TP_outlier"] = df["TP"] - 6.7
    df.loc[(8.3 < df["TP"]), "TP_outlier"] = df["TP"] - 8.3
    
    df["Alb_error"] = 0
    df["Alb_outlier"] = 0
    df.loc[(df["Alb"] < 4.4) | (5 < df["Alb"]), "Alb_error"] = 1
    df.loc[(df["Alb"] < 4.4), "Alb_outlier"] = df["Alb"] - 4.4
    df.loc[(5 < df["Alb"]), "Alb_outlier"] = df["Alb"] - 5
    
    # add tsne features
    header = ["latent{}".format(i) for i in range(2)]
    sc = StandardScaler()
    standarlized_data = sc.fit_transform(
        df.drop(columns=["disease"])
    )
    latent_features = pd.DataFrame(
        bhtsne.tsne(
            standarlized_data,
            dimensions=2,
            rand_seed=0
        ),
        columns=header
    )
    
    # add pca features
    header = ["pca{}".format(i) for i in range(3)]
    pca = PCA(n_components=3)
    pca_features = pd.DataFrame(
        pca.fit_transform(standarlized_data),
        columns=header
    )
    
    df = pd.concat([
        df, latent_features, pca_features
    ], axis=1)
    
    
    return df.iloc[:train_df.shape[0]], df.iloc[train_df.shape[0]:].drop(columns=["disease"]), ohe

学習

さて、学習に関してはハイパーパラメータをどう選ぶかも重要な要素になるのかと思いますが、私はそこについてはド素人です。
というか、そんなに興味が・・・ということで、Optunaに任せてみることにしました。

lightgbmのOptunaインテグレーションを使えば簡単にCVまでやってくれるのがわかったので、こちらを参考(丸パク)にしています。
ちなみに、うちの環境だと数回回すとJupyterKernelがクラッシュするという事象が発生していました。が、なんとか使えたのでそのままやってしまいました。

class TunerCVCheckpointCallback(object):

    def __init__(self):
        self.cv_boosters = {}

    @staticmethod
    def params_to_hash(params):
        params_hash = hash(frozenset(params.items()))
        return params_hash

    def get_trained_model(self, params):
        print(params)
        params_hash = self.params_to_hash(params)
        return self.cv_boosters[params_hash]

    def __call__(self, env):
        params_hash = self.params_to_hash(env.params)
        if params_hash not in self.cv_boosters:
            print(env.params)
            self.cv_boosters[params_hash] = env.model


def train(X, y):
    lgb_train = lgb.Dataset(X, y)

    # データセットの分割方法
    folds = StratifiedKFold(n_splits=4, shuffle=True, random_state=71)

    lgbm_params = {
        'objective': 'binary',
        'metric': 'auc',
        'verbosity': -1,
        'feature_pre_filter':False,
    }
    checkpoint_cb = TunerCVCheckpointCallback()
    callbacks = [checkpoint_cb]
    tuner_cv = lgb.LightGBMTunerCV(
        lgbm_params, lgb_train,
        num_boost_round=1000,
        early_stopping_rounds=100,
        verbose_eval=20,
        folds=folds,
        callbacks=callbacks
    )

    tuner_cv.run()

    print(f'Best score: {tuner_cv.best_score}')
    print('Best params:')
    
    p = tuner_cv.best_params
    # ※1
    p['num_iterations'] = 1000
    p['early_stopping_round'] = 100
    return checkpoint_cb.get_trained_model(p)

def predict(bst, x_test):
    return bst.predict(x_test, num_iteration=bst.best_iteration)

コード中に※1と記載しましたが、ここは参考にしたコードのままだと動かなかったのでこのように変更しています。
TunerCVCheckpointCallbackの__call__メソッド内にあるenv.paramsとtuner_cv.best_paramsで保存しているパラメータが異なるようで、そのまま渡すと該当するモデルが見つからずエラーになります。
もっと別の方法で、モデルの保存が必要かもしれません。

feature_importance

この結果のfeature_importanceを表示します。

pd.DataFrame(np.mean(bst.feature_importance(), axis=0), index=test_lgb_df.columns, columns=["importance"])

f:id:marufeuillex:20200905135916p:plain

少し考察してみます。

_errorは案外効いていないようです。この特徴量はその項目が正常範囲内なら0, それ以外なら1を取る変数でした。これは、おそらくデータから見て取れてそもそもの値から木に反映されやすいのだろう、と推測されます。

逆に**_outlier系はそこそこ効いているようです。これは正常値からの離れ具合を表しているので、遠ければ遠いほど異常になりやすいと考えられるので木としては使いやすかったのかもしれません。

pcaやtsneも結構効いてますね。

が、ここまで言ってきてなんですが、今回のコンペは何も考えずに最初からあった変数たちがかなり効いています。。
試してないですが、もしかしたらそのままの変数をいれてチューニング少しかければ0.92超えたのかも??なんてことも考えてしまいました。そうでないことを祈る。。

おわりに

今回は先日参加したBeginnerコンペについて書いてみました。
データがわかりやすいので脳筋な僕でもなんとか基準値を超えることができました。

これを弾みに、本格的なものも少しずつスコアを上げていけるよう頑張ってみます。