kaggle:PLAsTiCC Astronomical Classificationコンペのまとめ

# この記事の概要
2018/9/28～12/17で開催されたkaggleの「[kaggle PLAsTiCC Astronomical Classification](https://www.kaggle.com/c/PLAsTiCC-2018)」コンペに「Baba-horse」というチーム名で参加していました。結果は40th/1,100teamsとなんとか上位5%に入り所謂silver medalの称号を頂くことができました。
この記事では、どんなコンペだったか？の概要と自分たちの取ったアプローチ・上位者のアプローチとそこからの学びを、コンペに参加しなかった方もなんとなくわかる程度にまとめます。

![image.png](https://qiita-image-store.s3.amazonaws.com/0/123462/440ce860-dce4-c44e-bc94-12d84f8e074f.png)
# コンペの背景
※背景よりもコンペのことが知りたい！ という人はここは読み飛ばすことをお薦めします。
まず今回のコンペのhostをしている団体ですが、LSSTという大きな望遠鏡を作ってる&サイエンスを考えている人たちです。LSSTのWikiの内容をかいつまんで説明すると

> 3.2ギガピクセルの巨大なカメラを使用し、9.6平方度という超広視野を達成する。この広視野のおかげで、望遠鏡設置場所から観測できる全ての天域の写真をわずか3晩で撮り終えることができる。これを繰り返すことにより、これまで十分に観測できなかったような短い時間スケールで変動する天体 (超新星爆発や地球近傍小惑星など) を詳細に観測することが可能になる。

要するに**巨大な望遠鏡に超巨大なカメラを載せて超広い星空の写真をとって、大量の天体の写真を何度も取りましょう**、そして色々な発見をしましょう！という天文学の発展にとっては夢のようなプロジェクトのようです。

ただしまだこの望遠鏡で観測を開始しているわけではなく、2019～2020年くらいから観測が始まるようです。

なので、今回のコンペのデータは、このLSSTという大きな望遠鏡で取ったデータ…**ではなく【もしこの望遠鏡から観測データを取得したらこういうデータが取れます】というシミュレーションのデータ**になっています。

# コンペ概要
### 一言で言うと？→ 天体の多クラス分類です。
- そもそも「天体」と一口に言っても、様々な形態が存在しています。その中でも、今回のコンペでは「明るさが時間で変動する天体」の分類を行います。
    - 例えば、今回のコンペの例だと以下のような形態があるようです。
        - 変光星
        - 潮汐破壊
        - 超新星爆発
        - キロノヴァ
        - AGN
        - マイクロレンジング増光現象
        - 食連星
        - 太陽系外惑星
        - などなど…
    - 今回は15種類の形態の天体データを用意し、それぞれの天体がどの種類のクラスに属しているかを予測します
    - 所謂「多クラス分類」と呼ばれるタスクです
    
- 評価指標
    - multi-weighted-loglossで評価されます
    - ここでPijは天体iがクラスjに属する確率、Yijは実際に属していたら1,属していなかったら0を示し、Nは各クラスの天体数、Mはクラスの数、Wiは天体iのテストデータにおける重み（その天体が属するクラスの比率で重み付け）を示します。
    - testデータ内のクラスの比率の偏りが補正された指標になります。
    
    ![image.png](https://qiita-image-store.s3.amazonaws.com/0/123462/c61444c6-3151-f582-4278-da043034d6ab.png)


# 提供されたデータ
背景の部分で言及したように、LSSTという大きな望遠鏡で空を撮影した時のシミュレーションデータになります。ホスト側からは主に2つのテーブルデータが提供されました。

 1. 天体そのものに関するデータ
    - その天体が属するクラス（今回のtarget。後に言及する「class99」というクラスを除いた14クラスが存在しています。）
    - 天体の座標・地球からの距離（赤方偏移（photoz/specz）という単位で表現します）など
    - レコードは天体の数でユニークになっており、trainデータ:約7,800行×12列 testデータ:約350万行×11列となっています。
 2. 天体を観測したときの明るさに関するデータ
    - 各天体を観測したとき、時系列で明るさがどう変化したかのデータです。
    - 明るさ（flux）、明るさの観測誤差、観測時の光の波長（passband）などのデータで構成されています。
    - レコードは天体×passband×観測時間でユニークになっており、天体毎にだいたい100～300レコードほど存在します。よってデータ量は膨らみ、train:約140万×6列 test:約4.5億×6列のデータになります。

# このコンペで難しいと感じたこと
（注：私の個人的意見です。おそらくTOPの方々は「これはハードルではない」「むしろこれが難しかった」など全く別の視点などがあると思いますので、ご指摘・ご意見あれば適宜頂けると幸いです！）
### 1.ドメイン知識が無い人がテーマ・データを理解するハードルが高い
まず、最初のとっかかりとして、天文学に触れていない人はテーマの理解とデータの理解が一番の困難です。おそらく主催者サイドもそれを認識しており、そのギャップを埋めるために十分な天体の分類についての情報およびデータについての情報が提供されていました。その情報をしっかり読み込むことで（少なくともこのコンペに挑むにあたっての）ギャップはかなり解消されると思います。しかし、勿論そんなすぐに理解できる情報の量・質ではなかったと思います。おそらくここを深く読み込んだ参加者はそこまで多くなかったのではないでしょうか。

![shimekiri_report_hakui_man.png](https://qiita-image-store.s3.amazonaws.com/0/123462/894f2401-effe-0f9a-fc9c-b8007cb04b25.png)


### 2. trainデータの少なさとそこから起きる過学習
trainデータ（上述の1. 天体そのものに関するデータのtrain用データ）は8,000レコード弱で、そこから14クラスを学習します。クラスによっては極端に数が少ないクラスがあり、学習が難しいと感じました。何も考えずにパラメータチューニングするとどんどん誤差が小さくなるのですが、submitするととんでもなく過学習している、ということが発生していました。

![ai_kanabou_buki.png](https://qiita-image-store.s3.amazonaws.com/0/123462/a47fc16f-4bb2-f349-dbfa-d4451e44f80e.png)

### 3. testデータの多さから起きる計算リソースとの闘い
上述した2. 天体を観測したときの明るさに関するデータはとても量が多く、特にテストデータは約4.5億程度のレコードがありました。ここから複雑な特徴量の計算をすると時間がとてもかかるし、メモリもあっという間に食いつくします。ここの処理もおそらくコンペ参加者はそれぞれの工夫をした部分なのではないかと推察します。

### 4. テスト用データにだけ存在するクラス「class99」
学習用のデータに存在するクラスは14クラスですが、テスト用のデータにのみ存在する「class99」というクラスが存在していました。勿論データが無いので学習はできません。
後述しますが、ここのアプローチによってかなり精度に差が出ていたのかなぁと感じます。

# 自分たちのアプローチ
上記に挙げたそれぞれの課題にどうアプローチしたかを示します。
### 1.ドメイン知識→身も蓋もないですが、専門家とチームを組みました
ちょうど弊社の新卒の方で天文学を専攻している方がいて、その人とチームを組みました。
チーム構成は2人で、
- 専門的な知識をもとに特徴量とかのアイデアをくれるメンバー
- 基本的にコンペで使う特徴量を作ったり、モデリングとかを行うメンバー（私）
の構成で進めていきました。役割分担はかなりきちんとできており、うまく進行できたかなと思います。連絡は基本slackでやり取りしてました。

#### 使用した特徴量
上記のチームで、以下のような特徴量を作成しました。
    - バンド毎の明るさの最大/最小/平均など普通に思いつきそうなベーシックな特徴量
    - ピークの明るさから半分になる時間、1/6.3になる時間などドメイン知識がある程度必要そうな特徴量
    - また、コンペ終盤で、天体のphotozにあわせて各fluxが属するbandの補正を実施しました。
    - photozが大きい天体は波長が伸びるため、例えば観測上はband5でも実はband0として扱わなければいけないようです（by新卒astronomer。赤方偏移でググればわかるらしいですが、私はあきらめました。）。
    - この補正を実施して再度band毎の特徴量を算出しなおすとスコアを大きく上げることができ、ここでブロンズ圏からシルバー圏に押し上げてくれました。

### 2.trainデータの少なさと過学習 → binary_oofをした後に多値分類
trainデータをそのまま用いるのではなく、まずはクラス毎に0/1の二値分類で各クラスに所属する確率を予測した後、予測した確率を特徴量として改めて多クラス分類を実施しました。また、銀河系の天体と銀河系以外の天体はデータ上分離することが可能なのでデータを分けました。
所謂スタッキングの1種であり、自分のチームはこのアプローチでスコアを改善することができました。
また、そのまま多クラス分類をするとすぐに過学習するのでパラメータチューニングが難しいと思っていたのですが、二値分類のモデルに関してはチューニングで改善した部分がそのままスコアの改善にも反映されたので、二値分類の方はパラメータチューニングをゴリゴリ行い、多クラス分類のチューニングはpublic kernelと同じ設定でモデルを構築しました。

### 3. testデータのハンドリング→社内環境とGCPの利用
これも身も蓋もない話ですが、社内にメモリ64GBのサーバが複数台あり、勿論仕事上の使用が優先ですが空いている場合は仕事以外での用途で使用することが可能な環境が用意されています。基本的にはそこで作業していましたが、どうしてもメモリを食うような処理をしたい場合は別途個人的にGCPを立ち上げて作業していました。
ただし、基本的にはtestデータをまとめて処理するより分割して処理する方が当然効率が良いので、分割しながら社内サーバで作業していました。

### 4. class99は他クラスの確率で判断
これはカーネルで紹介されていたシンプルなアプローチで、結局最初から最後までこのアプローチを使用していました。以下の計算をベースとして、そこにクラスの重み補正を加えたものが共有されていました。

```math
\prod_{i=1}^{M}(1-pred_{i})
```
これだと例えば10クラス存在して各クラスで確率が0.1と予測されるとclass99は(0.9)^10で約35%となる一方、どこかのクラスの確率が1と予測されるとclass99は0%となります。

その他のアプローチとして、speczを加えたモデルやMLPのモデルを作ってensemble（単純な重み付平均）をしました。単体では精度が低いものの、重み付けして混ぜると精度向上に寄与してくれました。

これら1~4のアプローチで予測した結果、40th/1,100teamsとなりました。


ただし、私たちのチームのloglossが0.920だったのに比べ、TOP層は0.7~0.68とものすごいBig GAPが存在します。これらの方々はどういうアプローチを取っていたのでしょうか。


# TOP層のアプローチ
このコンペは上位者が取った手法がそれぞれかなり異なっており、色々なアプローチを取ることのできた良コンペであったと思います。上位の方が取っていたアプローチで特徴的だと感じたことを抜粋してまとめておきます。
### 1st Place Solution
https://www.kaggle.com/c/PLAsTiCC-2018/discussion/75033

この方は超新星爆発に関する研究者とのことで、ドメイン知識の豊富さで他を圧倒していた印象です。Discussionに書く内容も非常に専門的で、参加者にドメイン知識を提供してくれると共に自分でも圧倒的なスコアを残す、まさにスーパーマンのような存在でした。

- trainデータの拡張：testデータとtrainデータの明るさの傾向が違うので、testデータの傾向に合うようにtrainの天体1つにつき40パターンのデータを生成したとのこと
        - 具体的には明るさを変えたりphotozなどを変えたりしていたようです。何かの手法を使ったというより、ドメイン知識に基づくデータ生成という感じでしょうか。
   
- Gaussian Processによるライトカーブの予測：私はあまりGPに詳しくなかったのですが、以下の記事が理解の参考になりました。
[SVMより良いという噂のGaussian Processによる回帰を実装してみる](https://qiita.com/sergeant-wizard/items/d5a60d9f9ec08ef2d907)

- class99のLB probing：class99の確率をLBを見ながらチューニングした結果、class99の確率は特定のクラスの予測確率を重み付けした確率だということに気づいたとのこと。
- LightGBMのみで予測。アンサンブル一切無し（！！！）


### 2nd Place Solution
https://www.kaggle.com/c/PLAsTiCC-2018/discussion/75059

- 私の体感だと、参加者のメインモデルはLightGBMを使っている方が多い印象でしたが、この方はNeural Networkでメインモデルを構築していたようです
- 天体の明るさ（flux）とともにその誤差（flux_err）があることに目をつけ、flux_errの大きさに応じて特徴量として使うfluxを変えて数種類の特徴量データセットを作り、それを基にNNモデルを構築していたようです。
- また、class99のprovingやdata augmentation（ただし、1st placeのような補正は無し）やfluxの補正なども実施していたようです
- lightgbmとのアンサンブルで精度を上げたようですが、その際木の深さを浅くしたモデル（depth=2）を追加すると精度が上がったようです。

### 3rd Place Solution
3rd placeは日本人の方もいるチームになります。TOP3くらいでいつも争っており、（自分も参加しているのにも関わらず）最後はめっちゃ応援してました。
このチームは3名の方で構成されており、それぞれが異なるモデルを構築してからそれをマージしてスコアを上げています。

https://www.kaggle.com/c/PLAsTiCC-2018/discussion/75116
https://www.kaggle.com/c/PLAsTiCC-2018/discussion/75131
https://www.kaggle.com/c/PLAsTiCC-2018/discussion/75222


1. 1D-CNN
    - 様々な方法でサンプリングしたfluxの時系列(16channel×128 length)を用意
    - 最後の層で、別に作ったMLPの結果をconcatenated
    - MLPの結果は特徴量を変えて3種類作成
    
2. Catboost
    - データを銀河系と銀河系外に分けてモデルを構築
        - 銀河系モデル：2,000特徴量でモデル構築
        - 銀河系外モデル：1,000特徴量でモデル構築
        - 銀河系外モデル with specz：2,000特徴量でモデル構築
    - 上記の特徴量の組み合わせを変えて7×3パターンのモデルを作り、単純平均


3. Lightgbm
    - 使用している特徴量がとても専門的
        - sncosmoパッケージを利用したtemplate fitting（これが一番効いていた特徴量とのこと）
            - 私たちのチームもtemplate fittingには気づいていたのですが、気づいたのがコンペ終了1週間前で、計算リソース的にも完成するのがギリギリの時間なのでそれよりもモデルのensemble頑張ろう、という判断をしてしまいました…。
        - 時間（mjd）周りの特徴量
            - maxfluxのmjd - minflux（=detected:1）のmjdなど
        - 光度に関する変数
            - luminosity = (max(flux) - min(flux)) * distance ** 2 で算出
            - class90のデータのみを用いてpseudo labeling（データ増加）


また、これとは別に驚いたのはclass99 probingを徹底していたことでした。
最初のベースは私が上記で記載したclass99の確率を使っていたようなのですが、
あるクラスの確率の重みを変えるとスコアが上がることに気づき、そこからsubmitを重ねながら以下のformulaを導き出していたようです。すごすぎ…。
`extra_galactic_class_99 = 1.65 * (c42 + c52 + c62 + c95) ** 2.5 *  (1 - c95 - c62) ** 0.635'`
`galactic_class_99 = 0.2 * http://np.prod ([(1 - c) for c in classes])`


# （私の）学び
- ドメイン知識の重要さ
    - いかにデータの意味を知り、それに合わせた特徴量を作れるか。本当に基本のことなのですが、それがいかに大事かを痛いほど教わりました。
- 特徴量やデータの管理の重要さ
    - ここらへんはかなりなぁなぁで最後まで行ってしまいましたので、反省点です。「このsubmitでどんな特徴使ったっけ？」みたいなケースがたまに起きてしまいました（超反省）。特徴量やコード管理という基本のキをしっかり見直して次に繋げていきたいと思います！
- アイデアもそうだがそれを実装しきる覚悟やスピードが大事
    - template fittingやっぱやっとけばよかったです。ただし私たちには時間と覚悟とスピード感早く上記を抽出するスキルがおそらく足りませんでした。上位の方々の「やり切っている感」がすごいなと感じました。
    - 個人的には、これは1人で埋めるものではないかなと思っています。チームで補うことが大事なのでは、と感じるので次からは積極的にチーム組んでいきたいなと感じました。
- 特徴量の生成で、tsfreshは今後も汎用的に使えそう
    - コンペ途中で出てきたtsfreshのライブラリは、時系列データの特徴量をお手軽に出してくれるのでコンペでも仕事でも使いやすいし、非常に有用だなぁと思いました。

# 補足：社内の雰囲気
コンペ終盤は仕事を早めに終わらせてコンペに集中させてくれる環境があり、非常に助かりました。社内の有志メンバーで運用しているslackがあるのですが、そこにちょこちょこ進捗を共有すると応援が入るので、モチベーションもそこで維持することができました。
全く余談ですが、私の隣の席の後輩は「頑張ってください！」と応援しつつ、私がkaggleがっつりやってる隣で[**これ**](https://qiita.com/fujit33/items/dbfbd7a2aa3858067b6c)を作ってました。ある意味すごい会社です。



というわけでPLAsTiCCコンペについてまとめてみましたが、補足や訂正などあれば適宜ご指摘いただければ幸いです！