フリーランチ食べたい

Python/機械学習/データ解析/ソフトウェア開発などなど

Kaggleでよく使われるStacking/Blendingをheamy、Stacknetをpystacknetで高速に実装する

この記事は 機械学習工学 / MLSE Advent Calendar 2018 - Qiita の15日目です。

Stacking/Blendingは実装が面倒

  • Kaggleなどでよく使われるアンサンブル手法にStacking/Blendingがありますが自分で実装すると結構面倒です
  • またモデルの精度を測る際にK-Fold Validationを行うこともありますが、同じpredictを何度も発生してしまい、単純に行うと無駄な処理が多くなってしまいます。
  • 今回紹介するライブラリheamyでは、それらの問題点が解決されており、 抽象化されたAPIで簡単にStacking/Blendingを実装でき 、かつ、 内部でpredictされた結果がキャッシュされており、高速に計算を行うことができます
  • またheamyではWeighted Averageも簡単に実装できる仕組みがあります。
  • また、Stacking/Blendingに比べると浸透度は低いかもしれませんが、StacknetというStackingを複数階層に渡って行い、モデルでニューラルネットワークを作るような手法もあり、このStacknetもpystacknetというライブラリで簡単に実装することができます。
  • Stacking/Blendingに関しては下記のKaggleのドキュメントに詳しいです。

Kaggle Ensembling Guide | MLWave

  • 今回使ったコードは末尾に貼っています。

事前準備

データ

今回はscikit-learnに付属しているボストン市の住宅価格のデータセットを使おうと思います。

X, y = load_boston(return_X_y=True)
X_train, X_val, y_train, y_val = train_test_split(X, y, test_size=0.05, random_state=seed)

環境・精度計測方法について

環境は下記です。

$ sw_vers
ProductName:    Mac OS X
ProductVersion: 10.14
BuildVersion:   18A391

精度検証は10-Fold CVでmean_absolute_errorを使って行います。

ベースライン

単に使い方だけでなく精度や実行速度についても検証を行うため、まず単純な回帰モデルの実装と精度計測、速度計測を行います。 モデルは決定木、線形回帰、KNNと(決定木でもありますが)GBDTとしてCatBoostを使いました。

%%time
models_dic = {
    'random forest': RandomForestRegressor(n_estimators=50, random_state=seed),
    'linear regression': LinearRegression(normalize=True),
    'knn': KNeighborsRegressor(),
    'catboost': CatBoostRegressor(custom_metric=['MAE'], random_seed=seed, logging_level='Silent')
}

for name, model in models_dic.items():
    kfold = KFold(n_splits=10, random_state=seed)
    cv_results = cross_val_score(model, X_train, y_train, cv=kfold, scoring="neg_mean_absolute_error")
    print(f'{name} : {-np.mean(cv_results):.2f}')

精度

random forest : 2.26
linear regression : 3.38
knn : 4.29
catboost : 2.32

時間

CPU times: user 1min 47s, sys: 17.5 s, total: 2min 4s
Wall time: 49.5 s

ここではRandomForestが最も良い結果になりました。

heamy

それでは、heamyを使ってこれらのモデルをStackingしていきたいと思います。

github.com

必要APIのロード

from heamy.dataset import Dataset
from heamy.estimator import Regressor
from heamy.pipeline import ModelsPipeline

Stacking

%%time
# datasetを準備
dataset = Dataset(X_train, y_train, X_val) # X_testは今回使わないが入れないとエラーになる

# アンサンブルに使うモデルを定義
models = [
    Regressor(dataset=dataset, estimator=RandomForestRegressor, parameters={'n_estimators': 50, 'random_state': seed}, name='rf'),
    Regressor(dataset=dataset, estimator=LinearRegression, parameters={'normalize': True}, name='lr'),
    Regressor(dataset=dataset, estimator=KNeighborsRegressor, name='kr'),
    Regressor(dataset=dataset, estimator=CatBoostRegressor, parameters={'custom_metric': ['MAE'], 'random_seed': seed, 'logging_level': 'Silent'}, name='cr')
]

# pipelineを定義、2nd levelデータセットの作成
pipeline = ModelsPipeline(*models)
stack_ds = pipeline.stack(k=10, seed=seed)

# modelを作ってvalidation
stacker = Regressor(dataset=stack_ds, estimator=LinearRegression)
y_trues, y_preds = stacker.validate(k=10)

# 精度出力
cv_results = []
for y_true, y_pred in zip(y_trues, y_preds):
    cv_result = mean_absolute_error(y_true, y_pred)
    cv_results.append(cv_result)
print(f'stacking: {np.mean(cv_results):.2f}')

# X_testを使ってpredict
y_pred = stacker.predict()

精度

stacking: 2.14

時間

CPU times: user 1min 53s, sys: 18.4 s, total: 2min 12s
Wall time: 54.1 s

RandomForestよりも少しだけ良い結果になりました。 ポイントとしては計算時間もほとんど増えていないことです。 これは2nd levelデータセットの作成時にキャッシュが効いていることを意味しています。

Stackingの内部

上で見たようにheamyでは処理がかなりラップされており使い勝手は良いのですが、それぞれのオプション引数の意味が理解し辛いので少しだけ図解で解説します。

まず、Stackingでは2nd levelデータセットの作成を行います。ここでData Leakageを防ぐためにK-Fold CVを用いた手法を使います。 pipeline.stack(k=10, seed=seed) の部分ですね。図で表すと下のような感じです。

f:id:mergyi:20181215200828p:plain
2ne levelデータセットの作成

次に、結果のvalidationを行います。ここでは「predict&validate」になっていますが、純粋に予測したい場合は「predict」だけすることになります。 そしてこれを一般的なK-Fold Validationで精度検証しています。 stacker.validate(k=10) の部分ですね。

f:id:mergyi:20181215201220p:plain
2nd levelモデルの学習、結果のvalidation

このように、実際にやっている作業はかなり複雑なのですが、heamyではかなり単純なAPIで実装することができます。素晴らしい👏

Blending

Blendingもほとんど同じ実装で行うことができます。BlendingはStackingとかなり似ているのですが、2nd levelデータセットの作成でK-FoldでなくHoldout法を使うところが違いです。

%%time
# datasetを準備
dataset = Dataset(X_train, y_train, X_val) # X_testは今回使わないが入れないとエラーになる

# アンサンブルに使うモデルを定義
models = [
    Regressor(dataset=dataset, estimator=RandomForestRegressor, parameters={'n_estimators': 50, 'random_state': seed}, name='rf'),
    Regressor(dataset=dataset, estimator=LinearRegression, parameters={'normalize': True}, name='lr'),
    Regressor(dataset=dataset, estimator=KNeighborsRegressor, name='kr'),
    Regressor(dataset=dataset, estimator=CatBoostRegressor, parameters={'custom_metric': ['MAE'], 'random_seed': seed, 'logging_level': 'Silent'}, name='cr')
]

# pipelineを定義、2nd levelデータセットの作成
pipeline = ModelsPipeline(*models)
stack_ds = pipeline.blend(proportion=0.2, seed=seed)

# modelを作ってvalidation
stacker = Regressor(dataset=stack_ds, estimator=LinearRegression)
y_trues, y_preds = stacker.validate(k=10)

# 精度出力
cv_results = []
for y_true, y_pred in zip(y_trues, y_preds):
    cv_result = mean_absolute_error(y_true, y_pred)
    cv_results.append(cv_result)
print(f'blending: {np.mean(cv_results):.2f}')

# X_testを使ってpredict
y_pred = stacker.predict()

精度

blending: 3.04

時間

CPU times: user 11.2 s, sys: 1.83 s, total: 13 s
Wall time: 5.37 s

性能で見るとStackingより大分落ちてしまいました。 実際BlendingのStackingに対する特徴は実装がシンプルになる一方で過学習を起こしやすくなることなので納得のいく結果かもしれません。 ただ速度はかなり早いです。これは2nd levelのデータセットの作成でK-Foldを行っていない、かつ、ValidationのK-Foldでキャッシュが効いていることが原因な気がします。 もし他の要因や自分の計測方法に問題がありそうでしたら、教えてください。

Weighted Average

Weighted AverageはStacking/Blendingとは考え方が異なりますが、この手法に関しても簡単に行えるAPIが用意されているので紹介します。

%%time
# datasetを準備
dataset = Dataset(X_train, y_train, X_val) # X_testは今回使わないが入れないとエラーになる

# アンサンブルに使うモデルを定義
models = [
    Regressor(dataset=dataset, estimator=RandomForestRegressor, parameters={'n_estimators': 50, 'random_state': seed}, name='rf'),
    Regressor(dataset=dataset, estimator=LinearRegression, parameters={'normalize': True}, name='lr'),
    Regressor(dataset=dataset, estimator=KNeighborsRegressor, name='kr'),
    Regressor(dataset=dataset, estimator=CatBoostRegressor, parameters={'custom_metric': ['MAE'], 'random_seed': seed, 'logging_level': 'Silent'}, name='cr')
]

# pipelineを定義
pipeline = ModelsPipeline(*models)

# 最適な重みを探索
weights = pipeline.find_weights(mean_absolute_error)
pipeline_apply = pipeline.weight(weights)

# 精度を出力
cv_results = pipeline_apply.validate(scorer=mean_absolute_error, k=10)
print(f'weighted average: {np.mean(cv_results):.2f}')

# X_testを使ってpredict
y_pred = pipeline_apply.execute()

精度

weighted average: 2.16

時間

CPU times: user 2min 3s, sys: 21.2 s, total: 2min 24s
Wall time: 1min 4s

今回は精度としてStackingに負けてしまいましたが、RandomForestやCatBoost単体よりは良い結果になりました。 Weighted AverageがStackingよりも良い結果を出すことも多々あります。

pystacknet

次にStacknetを簡単に実装できるpystacknetです。

github.com

「そもそもStacknetとは?」という質問に関しては下のGithubリポジトリを見てもらうのが良いと思います。 元々Javaで実装されたライブラリでそれをPythonに移植したものがpystacknetになります。 簡単に言うと「Stackingをさらに積み重ねてニューラルネットワークのような形にしたもの」です。

github.com

実装

from pystacknet.pystacknet import StackNetRegressor
%%time
# アンサンブルに使うモデルを定義
models=[ 
    # 1st level
    [
        RandomForestRegressor(n_estimators=50, random_state=seed),
        LinearRegression(normalize=True),
        CatBoostRegressor(custom_metric=['MAE'], random_seed=seed, logging_level='Silent')     
    ],
    # 2nd level
    [
        Ridge(normalize=True),
        ExtraTreesRegressor(random_state=seed),
        XGBRegressor(random_state=seed)
    ],
    [
        LinearRegression(normalize=True)
    ]
]

# StackNetモデルを作成
model = StackNetRegressor(
    models, folds=10,
    restacking=False, use_retraining=True,
    random_state=seed, n_jobs=-1
)

# 精度出力
kfold = KFold(n_splits=10, random_state=seed)
cv_results = cross_val_score(model, X_train, y_train, cv=kfold, scoring="neg_mean_absolute_error")
print(stacknet: {np.mean(cv_results):.2f}')

精度

stackingnet :  2.20

時間

CPU times: user 30.3 s, sys: 18.8 s, total: 49.1 s
Wall time: 11min 31s

けっこう学習に時間がかかったのですが、精度はCatBoostと同じでした…残念。 ただコンペティションですごい成果を収めていることもあるモデルなのでチューニング次第ですね。 時間がかかったのはheamyと違いキャッシュが効いていないことも大きいと思います。 あくまで本家のStacknetをPython移植したものなので、そこまで完成されているわけではなさそう、という印象でした。

結果まとめ

今回はチューニングなどに力を入れていないのであくまで目安ですが、結果をグラフにまとめてみます。

f:id:mergyi:20181215211921p:plain
計測時間と精度の関係

こう見るとheamyを使ったStacking実装が有効そうな気がしますね。 もちろん色んなケースを試すべきだと思うので、もし全然異なる結果が出るケースなどあれば是非教えてください。

最後に

  • Stacking/Blending/Weighted Averageを簡単に、かつ、高速に実装できるheamyとStacknetを簡単に実装できるpystacknetを紹介しました。
  • ちゃんとStacking/Blendingを実装するとかなり面倒なので今後はどんどんheamy使っていこうと思いました。速度もかなり魅力ですね。
  • Stacknetは速度面ではラップトップで実際に活用するのは難しいかもしれません。

ちなみに

Stacking(Blendingも考え方は同じ)の概念に関しては以前記事を書いたのでこちらももし良ければ参照ください。 blog.ikedaosushi.com ※こちらは概念の理解を優先するために2nd levelのデータを作る際のCVを省略していますのでご注意

Gist

gist.github.com