🎃

[Kaggle] DSB2019: QWKをGBDTで直接最適化する

2020/11/09に公開

機械学習

kaggle

idea

1. 本記事について

私はKaggleのコンペ「2019 Data Science Bowl」（以下、DSB2019）において、Quadratic Weighted Kappaという少々扱いづらい評価指標に対して、GBDTのCustom Objectiveで直接の最適化を試みる工夫を行いました。

Private LB確定時のcheater排除の影響を受け、最終的にGold圏から外れてしまった経緯もあり、当初は公開する気をなくしてしまっていたのですが、以前からこういう技術記事を書いてみたいと思っていた中、突如やる気が起こったので勢いで書いてみることにしました。

ちなみに、DSB2019のデータセットは、どうやらホストの要請により削除されてしまっているようです。後々深く検証しようと思っていたので、それができなくなってしまったのは残念です。

最初に断っておくと、疑似データによる検証の結果、精度向上の面で本当に効果的なアプローチなのかを結論付けることは現状できていません。以下、少々長くなりますので、それをご承知の上で読み進めていただければ幸いです。

2. Quadratic Weighted Kappa

本コンペの評価指標は、順序付きの分類問題に用いられるQuadratic Weighted Kappa（以下、QWK）です。Evaluationのページによると、QWKの定義は以下の通りです。

\begin{aligned} w_{i,j} &= \frac{(i-j)^2}{(N-1)^2}, \\ \kappa &= 1 - \frac{\sum_{i,j} w_{i,j} O_{i,j}}{\sum_{i,j} w_{i,j} E_{i,j}} \end{aligned}

ここで、 $i$ は真値、 $j$ は予測値を表し、 $O_{i,j}$ は(真値, 予測値)の組が $(i, j)$ である度数です。また、 $N$ はターゲットの水準数です。 $E_{i,j}$ については後ほど説明します。

$\kappa$ は1以下の値を取り、1に近いほどよい指標です。

まず、 $w_{i,j}$ の分母 $(N-1)^2$ は、 $i, j$ に依存せず、 $\kappa$ の第二項の分母と分子の両方にかかりますので、無視できます。そうすると、 $w_{i,j}$ は単に真値と予測値の差の二乗、すなわち二乗誤差であるとみなせます。 $\kappa$ の第二項の分子は、それに $(i, j)$ の度数を乗じて総和しているので、残差平方和となります。

つまり、 $\kappa$ の第二項の分母さえなければ、 $\kappa$ の最大化は二乗誤差最小化と等価になるわけです。逆に言うと、この問題における曲者は、まさに $\kappa$ の第二項の分母というわけです。

$\kappa$ の第二項の分母についてですが、まずは $E_{i,j}$ について説明します。

$E_{i,j}$ は予測値の各クラスの割合をそのままに、真値と予測値の相関を完全に無くした場合に期待される $(i, j)$ の度数です。これは、予測値をシャッフルした際に期待される $(i, j)$ の度数と考えればわかりやすいでしょう。つまり、 $\kappa$ の第二項の分母は、予測値をシャッフルした際に期待される残差平方和であると解釈することができます。

ここまでのことから、（本来の意味はさておき）QWKという評価指標は、予測値をシャッフルした場合に期待される残差平方和に比べ、どの程度残差平方和を減少させることができているか、その比によって予測のよさを評価していると考えることができます。

ナイーブな例として、全てを同一のクラスで予測した場合は、予測値をシャッフルしたところで残差平方和は変わりませんので、 $\kappa$ の第二項は1となり、 $\kappa = 0$ となります。これは分類問題の評価指標の性質として直感に合っています。

では、意味のある予測をする場合は、どのように働くのでしょうか。

例えば、DSB2019のターゲットは0～3の4クラスでしたが、1や2に予測を集中させるよりも、0や3の予測を増やした方が、シャッフルした際に期待される二乗誤差は大きくなります。つまり、 $\kappa$ の第二項の分母は大きくなり、結果として $\kappa$ も大きくなる方向に働きます。

このことから、二乗誤差を最小化する場合に比べ、両端の値をより積極的に出力する方が、QWKにおいては最適になります。

実際本コンペで二乗誤差最小化で回帰を行った場合、真ん中辺りに予測値が集中するため、閾値を最適化して0～3に予測値を割り振る必要がありました。

3. 連続的な予測値に対するQWKの拡張

DSB2019のDiscussion等を見ていると、この問題に対するアプローチとして、

二乗誤差最小化の回帰問題として解く
ターゲットの順序は無視して分類問題として解く

のどちらかを採用した人が支配的だったように思います。

しかし私は、本コンペに取り組み始めて最初の数日間、どのように学習すればこの評価指標を最適化できるかを、文字通り紙とペンを使って考え続けました。具体的には、GBDTのCustom Objectiveで直接QWKを最適化するように学習することはできないか、ということです。

しかし、QWKは離散的な予測値に対する評価指標ですので、これをGBDTの上で最適化しようと思うと、まずは連続的な予測値に対しても適用できるような拡張を考えなければなりません。

QWKの定義式からは、一見すると $O_{i,j}$ や $E_{i,j}$ を計算することが必要に思えますが、実はそうでもありません。

既に述べたように、第2項の分子は残差平方和に等しいので、各予測値に対する総和として表せます。そして分母も同様に、各予測値に対して何らかの値を総和することで計算することができるのです。つまり、

\kappa = 1-\frac{\sum_{k=1}^{K} f(\hat{y}_k)}{\sum_{k=1}^{K} g(\hat{y}_k)}

の形に書けるということです。（ $i, j, N$ が既に使われてしまっているので、ここではデータに対する添字を $k$ としました。）

$f$ は予測値に対する二乗誤差ですが、 $g$ はどのようなものでしょうか？少し考えてみればわかりますが、これは予測値 $\hat{y}$ とテストデータにおけるターゲットのラベル割合によって決まる関数になります。

例えば、予測値を山勘で0としたときの二乗誤差の期待値は、

1^2 \cdot P(y=1) + 2^2 \cdot P(y=2) + 3^2 \cdot P(y=3)

で計算できることがわかるでしょう。

ターゲットのラベル割合は学習データから類推できます。DSB2019ではtrainとtestでラベルの分布が若干異なるため注意が必要でしたが、その点も考慮すると下図のような関数になりました。

ここまでくれば、これを連続的な予測値に対して拡張することは容易にできそうです。私は上図の4点を二次関数でフィッティングすることで分母を近似し、

\tilde{\kappa} = 1-\frac{\sum_{k=1}^{K} (\hat{y}_k - y_k)^2}{\sum_{k=1}^{K} ((\hat{y_k}-a)^2 + b)}, a \simeq 1.85, b \simeq 1.63

という関数によって、連続値を取る予測値に対してQWKを近似的に拡張しました。以下、これを拡張QWKと呼ぶことにします。

2023/2/14: 上式中の $a$ , $b$ は、それぞれターゲットの平均と分散に等しく、同時に近似ではないことがわかりました

4. GBDTによるQWKの最適化

QWKを連続値に対して拡張できたので、あとは第2項を抜き出して勾配を計算すれば、GBDTのCustom Objectiveに突っ込めそうです。少々複雑ですが、 $\tilde{\kappa}$ の第2項の勾配は以下のようになります。

\frac{2(\hat{y_k}-y_k)}{\sum_{k=1}^{K} ((\hat{y_k}-a)^2 + b)} - \frac{2 (\sum_{k=1}^{K} (\hat{y}_k - y_k)^2) (\hat{y_k}-a)}{(\sum_{k=1}^{K} ((\hat{y_k}-a)^2 + b))^2}

2023/2/11: 上式における第2項に誤りがあったのを修正しました

通常、XGBoostやLightGBMでは2階微分（hessian）まで用いるのですが、実は先ほど定義した拡張QWKは凸関数ではなく、2階微分は常に正とはならないため、hessianを使って最適化を行おうとするとおかしなことになってしまいます。

そこで私は、古典的なGradient Boostingと同様、2階微分は使わないことにしました。これは単に、全データに対してhessianを定数としてCustom Objectiveを定義するだけです。

私はRで本コンペに取り組んでいましたが、私が定義したCustom ObjectiveをPythonで書くと以下のようになります。

def qwk_obj(preds, dtrain):
    labels = dtrain.get_label()
    preds = preds.clip(0, 3)
    f = 1/2*np.sum((preds-labels)**2)
    g = 1/2*np.sum((preds-a)**2+b)
    df = preds - labels
    dg = preds - a
    grad = (df/g - f*dg/g**2)*len(labels)
    hess = np.ones(len(labels))
    return grad, hess

2023/2/11: 上記コードにおけるgradの計算式に誤りがあったのを修正しました

前述の数式と異なる点として、二乗誤差最小化とのアナロジーでf, gに1/2を付しているのと、データ数が多いほどgradの値が小さくなってしまうのを打ち消すために、gradにデータ数を掛けています。

また、gradを計算する前にラベルの両端で予測値をクリッピングしていますが、これは確信度が高い予測に対して、0を下回ったり3を上回ったりする予測値を許容するためです。

以上の工夫により、GBDTにおいてQWKを近似的に直接最適化することができるようになりました。この工夫は大きな差別化になると思い、本コンペを最後までやり切るモチベーションとなっていました。

5. 疑似データによる検証

本来はLate Submissionで検証を行いたかったところなのですが、データセットが削除されてしまったので、簡単な疑似データによる検証を行います。

疑似データは以下のPythonコードで生成しました。

import numpy as np
np.random.seed(0)
X_train = np.random(rand(100000), 10)
y_train = X_train.mean(axis=1)
y_train = np.floor(np.argsort(np.argsort(y_train)) / len(X_train) * 4)
X_train = X_train[:, :5]

これは、以下のような手続きです。

100,000のインスタンスに対して一様乱数を10個ずつ生成
インスタンス毎に10個の乱数の平均を計算
その平均が小さい順に0～3のクラスに4等分し、ターゲットとする
元の乱数10個の内、5個のみを説明変数として残す

シード値のみを変え同様の手続きでテストデータを作成し、目的関数が二乗誤差の場合と拡張QWKの場合で、QWKを評価指標として比較を行いました。

XGBoostとLightGBMでそれぞれ実験しましたが、簡単のためパラメータチューニングはせず、以下の設定で行いました。その他のパラメータはデフォルトです。

learning_rate: 0.1
max_depth: 5 (XGBoostのみ)
base_score: 1.5 (XGBoostのみ)
eary_stopping_rounds: 100

ここで注意いただきたいことがあります。色々と実験してわかったことなのですが、上記のbase_scoreの設定は今回の問題において重要です。なぜかというと、先ほども述べた通り拡張QWKは凸関数ではないため、初期値によっては上手く学習し切れず、精度を落とすことがあるということがわかったからです。そのため、予測値の初期値は真ん中当たりに指定した方が無難です。

LightGBMには、（私の知る限り）base_scoreを設定する方法が無く、予測値の初期値は0からスタートしてしまいます。（通常の回帰や分類の場合はターゲットの平均値を初期値としてくれるのですが、今回のようにCustom Objectiveを使う場合、なぜかそうしてくれないのです。）

そこで、今回の実験でLightGBMに拡張QWKを適用する際は、一旦ターゲットの範囲を $[-1.5, 1.5]$ にシフトしてから学習し、最後にシフトし直すことで評価を行いました。

以上色々と述べましたが、肝心の実験結果は以下の通りです。

Method	Objective	test QWK (*1)	test QWK (*2)	分割閾値
XGBoost	二乗誤差	0.5806	0.6335	[0.928, 1.516, 2.098]
XGBoost	拡張QWK	0.6350	0.6351	[0.462, 1.492, 2.506]
LightGBM	二乗誤差	0.5797	0.6348	[0.896, 1.520, 2.112]
LightGBM	拡張QWK	0.6352	0.6347	[0.524, 1.480, 2.496]
*1: 予測値を単に四捨五入した場合
*2: 交差検証（5-fold）で求めた最適な分割閾値を適用した場合

パラメータチューニングを行っていないため、XGBoostとLightGBMとの間の比較にはあまり意味が無いことにご注意ください。

まず、XGBoostとLightGBMの両方とも、二乗誤差の場合は閾値を最適化しないと精度が出ないのに対し、拡張QWKの場合は単に四捨五入でクラスを割り当てた場合とほとんど差異が見られないことがわかります。期待通り、拡張QWKによってQWKを直接最適化できていると言えそうです。

では、閾値を最適化した上での精度はどうでしょうか。XGBoostは拡張QWKによりQWKが幾分改善したのに対し、LightGBMでは残念ながら向上が見られませんでした。

DSB2019ではXGBoostを使っていましたが、ローカルで検証した限りでは、拡張QWKによって明らかに精度向上が見られました。ただ、LightGBMだとどうなるかは検証できていませんし、ハイパーパラメータに依存する部分もあるかもしれません。

今回は単純なデータセットだったからかもしれませんが、これが一般的に有効なアプローチと言えるかどうかは、もう少し検証が必要そうです。

2023/2/11: PlayGroundコンペにおける検証コードを公開しました
https://www.kaggle.com/code/rsakata/optimize-qwk-by-lgb

6. 最後に

評価指標に合わせて学習のやり方を調整するというのは、あまり本質的ではないとする見方もあるかもしれません。

しかし、機械学習モデルがビジネスにもたらすインパクトを評価する際、そのインパクトの大きさが数学的に扱いやすいRMSEなどの指標に帰着するというケースは、むしろ多くはないと思います。

例えば、売上予測のタスクにおいて、誤差金額の絶対値がビジネスインパクトとして現れるのであれば、RMSEよりもMAEを最適化するような工夫をするべきでしょう。

QWKが現実的に意味がある指標かどうかはさておき、評価指標に合わせて学習のやり方を工夫するというのは、決して無意味なことではないと思います。

以上、長くなってしまいましたが、何かの参考になれば幸いです。