nykergoto’s blog

機械学習とpythonをメインに

Kaggle Days Tokyo のオンサイトコンペに参加しました! #kaggledaystokyo

Kaggle Days Tokyo で開催されたオンサイトコンペに参加してきました!! 結果としては全体 88 チーム中 private で 56 位という悔しさの残る結果になりました。が同時に反省点と学びもとても多い素晴らしいコンペだったので、感想兼反省文を書いていこうと思います。

f:id:dette:20191214155803p:plain

Kaggle Days Tokyo は2日間ありましたが、僕は体調不良のため2日目しか参加できませんでしたので1日目のことに関しては他の方の記事を参照していただければと思います。

どんなコンペだったか

日経電子版のログデータをもとにして、閲覧しているユーザーの年齢を当てるというタスクです。メインのテーブルには

  • ユーザーid
  • 見ているデバイス情報
  • 記事のid

などがありこの記事 id が記事データと紐づくような構造になっています。記事データには

  • タイトル (1/2/3)
  • キーワード
  • 本文
  • 記事の発行された日時

などがありました。

基本的な戦略

僕は弊社でインターンをしてくれている もーぐりくん と一緒にチームで参加しました。 僕はオレオレのフレームワーク vivid があり初速がある程度出せるだろうという見積もりがあったので、僕が基本的なモデリングを行い feature importance など特徴量の重要な話などは共有しつつ、最後にマージしましょうという作戦を取りました。

NYK510 TimeLine

ML Bearさん が当日の timeline を書かれていて、後で反省するのに良いと思ったので僕も真似します。 正直かなり切羽詰まっていたのでほぼ覚えていないのですが git の message とともに振り返っていこうと思います。

  • 10:35 initial commit
    • atma-cup #2 のリポジトリをコピペしてスタートしました。
  • 11:22 [update] version-1 feature
    • データの構造と submission すべき情報を見て submission までの雛形を作成していました。このときにつかっているのはメインのテーブルだけでした。
  • 11:50 [update] first submit
    • 最初のサブミット. 単体モデルでは Objective="poisson" の LightGBM が一番良かったためそれで submit しました. この日で一番余裕があったタイミングだったと思います。
    • 一旦 10 モデルぐらい作成して ridge 回帰も作っていたのですがそれぞれのモデルのチューニングが適当すぎたのか CV と LB の差が激しく撃沈
  • 12:01: FastText のモデル作成完了
  • 13:10: [update] add pseudo labeling
    • pseudo labeling をやってみたかったので実装してみていました。盛大にバグって精度が出ませんでした。このあたりから焦り始めます。
  • ~ 14:30 [fix] bug
    • NN を無理やり入れようと苦戦して bug と戦う羽目になりました。また optuna での tuning を回しながら特徴を作ろう、であったり記事のテキストを整形しよう、みたいな欲張りをしてそちらでも bug と戦っていました。この間の進捗は虚無でした。
    • このあたりで一回もーぐりくんとMTGをして正気を若干取り戻し、記事情報を少し入れました。またもーぐりくんアイディアの特徴も入れてスコアは良くなりました。
  • 15:17: [update] swem enbedding
    • FastText で作った特徴量で記事情報の SWEM をしました。この時点でかなり上位陣と差をつけられていて差分がわからず混乱に陥っていました。(後で判明しましたがこのときロジックのミスでちゃんと本文情報を使えていなかったみたいです)
  • 16:47 [update] genre
    • 記事のジャンル情報を入れてみて若干のスコア改善。もすぐに抜かされまくるという状況。
  • ~ 18:30 気づくとコンペはおわっていた

ふりかえり

悪かったこと

チームのメリットを活かせなかった

最初はまだしも、途中からは自分のことでいっぱいになり、全くチームとして機能していませんでした。

基本的に今回僕は独創的な特徴量は考えられませんでしたが、すこしの議論しかできませんでした。 一方でもーぐりくんはユーザーがどういう気持で閲覧するかなどを考慮した集計方法などを提案してくれていて、非常に頼もしかったです。 そのアイディアのパワーを生かしきれなかったのは反省です。

実装の正確性

今回は家のPCに ssh 接続して分析を行っていたのですが noetbook のバグ? かなにかで jupyter 上のブラウザからささっとファイル閲覧ができなくなり、最終的な output の特徴量をエクセルで確認するような作業をふっとばしていました。 これによって、よくある「概ねはあっているが、1行違うために生成される特徴量がおかしい」というパターンに陥って、性能がでないという自体になっていました。

次の日コードを見て気づいて修正すると Private で 11.70366 → 11.41330 (13位相当) になりました。悔しい。

知っている内容を使えなかったこと

今回の solution ではほとんどすべてのチームが TargetEncoding を使ってカテゴリの埋め込みを行っていて、どのチームもとても良く効いたとのことでした。 僕はリークを恐れて最後まで「TargetEncodingを使う」ことが頭の中の選択の一つにも上がっておらず、ここで大きな差をつけられてしまった感があります。

ではこの情報に僕が到達できなかったかというと否で、なんなら atmaCup#2 でも target encoding の使い方が肝でしたし分析コンペLT会のハクビシンさんの発表も、なんなら前日の Jack さんの発表でも TargetEncoding 大事だよーということは言われていて twitter などでも耳にしていましたから「TargetEncodingの実力を過小評価していた」自分の落ち度です。

hakubishin さんのスライド: Target Encoding はなぜ有効なのか

speakerdeck.com

Jack さんのスライド: How to encode categorical features for GBDT

speakerdeck.com

ちなみにメインテーブルのカテゴリ変数に対して TargetEncoding を投入すると Private で 11.70366 → 11.25799 (3位相当) になりました。悔しいなあ。

良かったこと

もーぐりくんがチームメイトだったこと

多分僕が二人のチームだったらもっと悲惨でした。ありがとう。

vivid で完走できたこと

今回のコンペで与えられていたテーブルは予測する際に一度ユーザーで集計をする必要があるデータでした。 そのようなデータは vivid を作ったときには全く想定していなかったのですが特徴変換の部分とそれを pipeline 的に行う部分を分離して実装していたおかげで対応できたので構造化がある程度できていたのかなと思っています。 またある程度の初速を出せたのは(あとでめちゃ抜かされてはいますが)うりではあるのでそれも良かった点でしょうか。

使っているうちに出てきたダメポイントについても、体を張ったテストだと思って、全部 issue にしてより強くしていきたいと思います。

たくさんの Kaggler と一緒にオンサイトコンペに参加できたこと

これは何者にも変えられない体験でした。特に夕方のある程度形勢が定まってきたような段階でも皆黙々と作業に打ち込んでいる様子は流石だなあと思っていました。 同じ場に参加できたこと、大変嬉しく思います。

またコンペのタスク・データについても非常に噛みごたえのある面白いデータでした。提供いただいた日経さん、また設計を主に担当されたであろう u++ さんに感謝です。

学び

精神的なこと

いつもは開催者側にいるのでなんとなーく大変なことはわかっていましたが、時間制限付きのコンペの大変さは僕の想像を遥かに上回っていました。*1

もとから知っていはいましたが僕がとても予想外・いつもと違う状況に弱いことを改めて知ることができました(悲しいことにほとんど覚えていないのですが夕方からコンペ終了までの僕はかなり狼狽していたと思います)。 そもそも、そういう状況になることを思って行動していない & 腹をくくって思考を切り替えられないのは分析コンペだけではなく必要な能力だと思います。正直ちょっとどうやって鍛えたら良いのかがわかっていませんが、精進します。

追試験

勉強のため上記で触れたような WordEmbedding・Pseudo Labeling*2 まわりなどのバグの修正 & TargetEncoding の追加等の修正など行い LateSubmission をしてみました。 マシンパワーの関係で LightGBM の SingleModel しか試せていませんが 一番良かったのは Poisson の LightGBM で Private 11.17052 (2位相当) となりました(SeedAveraging や Stacking など行えばまだ上がるかもしれないです)。

五月雨ですがやったことで効いた・効かなかったはざっくり以下のような感じになりました。

  • 効いたこと
    • target encoding
    • pseudo labeling
      • やるごとにじわっと良くなる (private で 0.02程度?)イメージでした。3回ぐらいまでは有効でした。(それ以上はサチって止まってしまう感じ)
    • SWEM (simple word embedding) での記事タイトルや本文、タイトル、キーワードの埋め込み
    • Objective=RMSE 以外で解く
      • 自分の範囲内では Poisson が一番良かったです
      • いま Log 変換を試していないことに気が付きました。後でやる。
    • ユーザーと記事との関係から記事の embedding の作成
    • seed averaging
    • キーワード全体での count で置き換えて集計 (mean/sum/std)
    • ユーザーごとのアクセス時間のヒストグラム特徴量 (もーぐりくん案)
  • (自分の範囲内では)効かなかったこと
    • one-hot 化して多クラス分類として解いた後に期待値に直す
      • 画像からユーザーの年齢を当てるタスクの論文で上記の方法のほうが regression よりも有効だった旨が書いてあったのを思ってやってみましたが余り効果ありませんでした。*3
      • そもそも学習させたモデルが LightGBM だったので駄目だったのかも知れません
  • 効果が微妙だったこと
    • カテゴリ変数の one-hot encoding
      • 水準数がとても多いカテゴリも存在していたので Target Encoding のほうが効率もよくあえて one-hot する必要はなかったかなと思っています。

コードは以下においてありますのでもしよろしければ。vivid 使ってます 😌

github.com

やっていて、当日の8時間の中でこれをやってしまう上位のチームのパワーに圧倒されるばかりでした。しかも優勝チームの pocket さんに話を伺ったとき、全部スクラッチで書いてーという話をされていて、圧倒的な自力の差を感じました。

最後に

このような機会を頂いたKaggle Days Tokyoの運営の方々, データ提供頂いた日経様, スポンサード頂いている会社様, ありがとうございました。 次回以降も機会あれば是非参加したいです!

参考文献

naotaka1128.hatenadiary.jp

py2k4.hatenablog.com

*1:あらためてatmaCupもそうですがオンサイトで勝つ方のタフさを感じました

*2:限定的ですが Pseudo Labeling は効果がありました

*3:ちょっと読んだのが昔なのでうろ覚えですがたぶん DAGER: Deep Age, Gender and Emotion Recognition Using Convolutional Neural Network この論文

RMSE を Fold ごとに取ると全体の値より小さくなる証明

この記事を書く前に twitter でお話をしている流れで、まますさんに的確な証明を頂くことができました! 証明にはこちら RMSE.pdf - Google ドライブ からアクセスできます。(まますさんありがとうございましたmm)

そもそも

この記事のお題は RMSE を Fold ごとに取ると全体の値より小さくなる証明をやります ということです。

これをやろうと思ったきっかけは #かぐるーど での kaggle本の本読みです。 前回は第5章だったのですが、その5.2.2で次のような記述があります。

クロスバリデーションでモデルの汎化性能を評価する際は、通常は各foldにおけるスコアを平均して行いますが、それぞれのfoldの目的変数と予測値を集めてデータ全体で計算する方法もあります。なお、評価指標によっては各foldのスコアの平均と、データ全体で目的変数と予測値から計算したスコアが一致しません。例えば、MAEやloglossではそれらが一致しますが、RMSEでは各foldのスコアの平均はデータ全体で計 算するより低くなります。

要するに K-Fold それぞれの rmse 平均値と、データセット全体での rmse の値だとデータセット全体のほうが大きい (K Fold のほうが良いように見積もられてしまう) という話です。たしかに Root をとる操作を毎回やるのと、全体で合わせた後やるのだと前者のほうが小さい値になような感じはしますよね。

これ一般的に示せるかなーという議論があり、僕が「関数の凸性とイェンゼンの不等式でいけますよ」と言ったところじゃあやってほしい!と言われたのが当エントリの経緯になります。

せっかくなので簡単にですが凸関数とイェンゼンの不等式にも触れつつ、お話できればと思っています。

NOTE: 若干細かい定義域についてやイェンゼンの不等式の導出についてなどは省略していますので、それらは別途文献など見ていただけば幸いです。

凸関数とは

下準備として、凸関数というのを定義します。凸関数というのは色々な定義がありますが、以下を満たすような関数 $f: \mathbb{R} \to (-\infty,+\infty]$ のことです

$$ f(t x_1 + (1 - t)x_2) \leq t f(x_1) + (1-t) f(x_2) $$

ただし $x_1, x_2$ は任意の実数 $\mathbb{R}$ の点で $t$ には 0以上1以下の制約がついています。

f:id:dette:20191206013409p:plain
wikipedia 凸関数より引用

要するに $x_1$ と $x_2$ の内分点での $f$ の値と最初に $f$ で計算してしまってから $f(x_1)$ と $f(x_2)$ の内分を取るのとだと、後者のほうが大きいような関数、って言うことです。

また $f$ が微分可能な場合 $f$ の二階微分 $f'' \geq 0$ であることと上記の等式は同値になります。

イェンゼンの不等式

これをちょっと発展させて内分点の部分を2つ以上の点に拡張したのがイェンゼンの不等式です。

イェンゼンの不等式は上記の式と同じく特定の関数 $f$ が凸関数である必要十分条件を表した式で, $f$が凸ならば任意の自然数 $n$ と$\sum_{i=1}^n p_i = 1, p_i \geq 0$ を満たすような $p_i$ に対して次の式

$$ f(p_1 x_1 + p_2 x_2 \cdots + p_n x_n) \leq p_1 f(x_1) + p_2 f(x_2) \cdots + p_n f(x_n) $$

がなりたつ、という定理です。

たとえば$n=2$の時を考えてもらうと先ほどの凸関数の定義そのままであることはすぐわかると思いますので、凸関数の定義を変数 $n$ 個の場合に拡張したようなイメージです。

RMSE を考える

RMSE とは入力とラベルの誤差の2乗和を $M$ とした時に

$$ {\rm RMSE}(M) = M^{\frac{1}{2}} $$

で計算される値です。これの二階微分を考えると

$$ {\rm RMSE}''(x) = - \frac{1}{4} M^{- \frac{3}{2}} < 0 $$

です。すなわちRMSEの二階微分は常に負の値となります。これは凸関数と全く正反対の性質で一般に凹関数 (concave) と呼ばれ先のイェンゼンの不等式とちょうど不等号が反対の不等式が成立します。

K-Foldしたときの RMSE

今Fold を$ K$ 個に分割して、それぞれが $n_k$ 個のデータを持っているとします。(データセット全体では $N$ 個とします。) この時各 Fold での MSE (Mean Squared Error) を $M_k$ とすると Fold ごとのデータの数で重みづけた ${\rm RMSE}_{\rm fold}$ は

$$ {\rm RMSE}_{\rm fold} = \sum_{k=1}^K \frac{n_k}{N} \sqrt{M_k} $$

となります。一方で通常の RMSE に関しては

$$ {\rm RMSE} = \sqrt{\frac{1}{N} \sum_{k=1}^K n_k M_k} = \sqrt{\sum_{k=1}^K \frac{n_k}{N} M_k} $$

となります。ここで $M_k$ に $n_k$ をかけているのは $M_k$ が既に Mean Squared Error なので要素の数を掛けて和にになおして全体の $N$ で割算をするためです。

ここで $p_k = n_k / N$, $f(x) = \sqrt{x}$ と考えると $f$ は凹関数でかつ $\sum p_k = 1$ ですのでイェンゼンの不等式が用いることが出来て

$$ {\rm RMSE} = \sqrt{\sum_{k=1}^K \frac{n_k}{N} M_k} \geq \sum_{k=1}^K \frac{n_k}{N} \sqrt{M_k} = {\rm RMSE}_{\rm fold} $$

が成立します。即ち fold ごとで RMSE を計算して重み付きの平均を取った値のほうが、データセット全体での RMSE の値より小さくなることがわかりました。

RMSE 以外でも…

上記の証明を追っていただくと分かるようにこの証明はロス関数の値がデータごとに計算できること、及びそれをデータセット全体の平均したあとに凹関数に代入する、という構造が保たれている限り同様の議論をすることが可能です。

ですので Log を取ってから Root を取る RMSLE (Root Mean Squared Log Error) なども同様の議論が可能です。

参考文献

以下は本記事を書くにあたって使用した凸関数に関する話題や凸最適化に関する日本語の参考文献です。

この記事は kaggle その2 advent calendar 2019 の記事です。

分析コンペLT会でLTをさせてもらいました!!

Kaggle 旅行

2019/11/30 に行われた分析コンペLT会にLT枠として参加させていただきました。😄

kaggle-friends.connpass.com

僕は普段大阪で仕事をしているのでもともと発表はおろか参加する予定はなかったのですが(行きたいとはめっちゃおもっていた)、いつもかぐるーどでお世話になっているカレーさんに「発表枠あけて待ってます！！」(原文ママ)と嬉しいオファーを頂いき、ちょっと東京旅行兼ねて発表枠として参加しました。

内容に関しては俵さんが素敵にまとめて頂いているのでそちらを参考にしていただければと思います。

tawara.hatenablog.com

発表内容

僕の発表は「初手が爆速になるフレームワークを作ってコンペ設計した話」というタイトルで、自分が作っている https://gitlab.com/nyker510/vivid という機械学習用のフレームワークとそれを使ってコンペ設計が楽になったよ、ということを話しました。

speakerdeck.com

Vivid について

Kaggle でもそうですがお仕事でも当然特徴量とモデルのバージョンを管理するのはとても大切です。なんですが僕はとても大雑把な人間なので、良く「このファイルどのスクリプト or Notebook から出てきたんだっけ…」ということに陥っていました。

そこで「動かすだけで勝手にログもモデルもバージョンも保存してくれるやつがあったらいいなーというかないと僕は無理だな」と思い、自分で色々と試行錯誤をして出来上がったフレームワークです。

大きなコンセプトというか特徴は以下のような感じ。

  • 必要なことだけを書くので良い
    • 基本は勝手にやってくれる
    • k-fold の split をして oof を計算するコード、みたいな定形処理は全部 vivid にお願いして、プロジェクト固有のコードに集中できる用に
    • ログの出力やモデルの保存なんかも勝手にやってくれる用に
  • テンプレート的な特徴作成の提供
    • 毎回 count encoding のコードを書くのは良くないのでそれもやってもらう
    • 特徴量をある粒度 (atom とよんでいます) とその集合体 molecule で管理して versioning する機能とかもあります
  • スタッキング・アンサンブル対応
    • 対応というよりは、一気通貫に出来るというのが売りです。
      (全部一気につながって作るので、前の run で作った特徴で学習していて予測するとバージョン違いで精度が出ない・カラムの数が違う、といった悲しいミスを防ぐ)

立ち位置としては sickit-learn よりも更に盛ったフレームワークという感じでしょうか。(webのDjango的な)

オレオレフレームワークに過ぎないのでコードもアレだしドキュメントもないしで正直めちゃ恥ずかしいんですが「こうやったらいいんじゃないか」とか「僕はこうやってますー」みたいな知見がでてお話できればとてもうれしいです。;)

ソースコードはこちら https://gitlab.com/nyker510/vivid からアクセスできます。 pip なら

pip install git+https://gitlab.com/nyker510/vivid

でインストールできますのでちょっと触っていただけるとめちゃ喜びます

atmaCupについて

スライドで atmaCup について触れたのですが結構な方に知っていただいてとても嬉しかったです! 次回第3回も鋭意企画中ですのでぜひご参加いただけるとこちらもとてもとても喜びます😆

感想

やはり発表すると思っていることがまとまるので、発表者が得るものが大きいなととても感じました。今回もこの発表のために vivid をちょっと直して(そしてバグを見つけ😌)、スライドにするために自分の考えを一度整理して、という作業をおこなったのでいま自分がやっていることの意味合いがクリアになって非常に良い経験でした。

また発表のことについて質問してディスカッションできたり、僕も他の人の発表について聞けて(そしてどれもめちゃくちゃクオリティが高い!!) 質問できたりというのは実際にあって喋る良い点だなあと思いました。

今後も機会があれば是非参加したいと思います! 企画頂いたかぐるーどの皆様、会場提供頂いた日経BP様、関係者の方々本当にありがとうございましたmm

f:id:dette:20191203003118j:plain
次の日はぶらぶら観光していました(写真はLT会の次の日に言った三菱一号館美術館です)。仕事でない東京もいいものですね。

PRML の本読みをしています @section3

PRML Python 機械学習 numpy

最近(というか今日から) 会社でPRML勉強会をやっています。ふつう第1章からやるのが普通ですがPRMLはちょっと重たいので息切れすると良くないよねということでいきなり第3章から始めるという方針をとっています。*1

今回は僕が担当で、主に

  • 線形モデルの導入
  • 正則化 (なんで Lasso はスパースになるのかの話)
  • バイアス・バリアンス

の話をしていました。たぶんこの話の中で一番一般的に使われるのは「バイアスバリアンス」の話だと思っていて、改めてまとめてみて個人的に良かったです。ざっと今日話した内容含めてまとめると以下のような感じかな?と思ってます。

f:id:dette:20191112005503p:plain

前提

  • $D$: データの集合
  • $E_D$ データ集合での期待値を取ったもの
  • $y$: 予測関数
  • $E_D[ y(x;D)]$ データの集合で期待値を取った予測関数。いろんなデータでモデルを作ってアンサンブルしてるイメージ
  • $h$: 理想的な予測値

バイアス (モデルの傾向)

  • データを沢山取ってアンサンブル平均した $E_D[ y(x;D)]$ が理想的な関数 $h$ にどれだけ近いか
  • ロジック自体の傾向が強い(biasが強い)と大きくなる
    • [低] 表現力が大きく多様性があるモデル
    • [高] 表現力が小さいモデル (いくら平均しても理想に近づけない)
      • たとえば常に 1 を返すようなモデル $y(x) = 1$ などはハイバイアスです。
      • 癖が強い(データに合わせる気がない)というイメージ

バリアンス (Variance: モデルの分散)

  • 今回のデータでの予測値と平均的なデータでの予測値の差分
  • 毎回のデータでどれぐらいばらつくか
    • [低] データによって予測が変わらないもの
      • たとえば学習を全くしない y=1 を常に返すアルゴリズムは Variance=0 です。反対にバイアスはめっちゃ大きい。
    • [高] データに依存して予測値が良く変わるモデル
      • たとえばパラメータ過多の最尤推定のようにデータに過剰に fitting してしまうアルゴリズムはハイバリアンスです

おまけ

バイアスバリアンスの説明に使われている Figure 3.5 にあたる絵を書く python script を書いてみました。

自分で書くと色々と考えるので楽しいですね。

import numpy as np

import matplotlib.pyplot as plt

def gaussian_kernel(x, basis=None):
    if basis is None:
        basis = np.linspace(-1.2, 1.2, 101)
    
    # parameter is my choice >_<
    phi = np.exp(- (x.reshape(-1, 1) - basis) ** 2 * 250)
    
    # add bias basis
    phi = np.hstack([phi, np.ones_like(phi[:, 0]).reshape(-1, 1)])    
    return phi

def estimate_ml_weight(x, t, lam, xx):
    basis = np.linspace(0, 1, 24)
    phi = gaussian_kernel(x, basis=basis)
    w_ml = np.linalg.inv(phi.T.dot(phi) + lam * np.eye(len(basis) + 1)).dot(phi.T).dot(t) # bias があるので +1 しています
    xx_phi = gaussian_kernel(xx, basis=basis)
    pred = xx_phi.dot(w_ml)
    return pred

n_samples = 100

fig, axes = plt.subplots(ncols=2, nrows=3, figsize=(10, 12), sharey=True, sharex=True)

for i, l in enumerate([2.6, -.31, -2.4]):
    ax = axes[i]
    preds = []
    for n in range(n_samples):
        x = np.random.uniform(0, 1, 40)
        xx = np.linspace(0, 1, 101)
        t = np.sin(x * 2 * np.pi) + .2 * np.random.normal(size=len(x))
        pred = estimate_ml_weight(x, t, lam=np.exp(l), xx=xx)
        
        if n < 20:
            ax[0].plot(xx, pred, c='black', alpha=.8, linewidth=1)

        preds.append(pred)

    ax[1].plot(xx, np.sin(2 * xx * np.pi), c='black', label=f'Lambda = {l}')
    ax[1].plot(xx, np.mean(preds, axis=0), '--', c='black')
    ax[1].legend()
    
fig.tight_layout()
fig.savefig('bias_variance.png', dpi=120)

f:id:dette:20191112002247p:plain
バイアス・バリアンスのやつ

*1:若干荒業感がありますが、一人以外は一度は読んだことがあるというのとやはり息切れが怖いのでちょっと先に応用が出てきそうなところから取り組んでいます