NFL 1st solutionの再現にチャレンジしました！

Kaggle Days Tokyoお疲れ様でした！これだけたくさんのKagglerが集まる機会はなかなかないので、とても楽しい2日間でした。

ということで、Kaggle Daysでもチームを組んだThe Zooのアメフトコンペのソリューションを再現すべく、これまでやった内容を投稿します。

https://www.kaggle.com/c/nfl-big-data-bowl-2020/discussion/119400#latest-689317

アメフトコンペPublic 1st Solutionの再現チャレンジ

モチベーションとしては、アメフトだけでなくサッカーやバスケなど、オフェンス、ディフェンスが入り組んで、プレイヤーの位置も常に動き続けるようなスポーツの解析全般に応用できると考えたからです。

結論から言うと、完璧な再現ができていませんが、せっかくなのでやった内容を順に説明します。

コンペ概要

コンペページはこちらです。

タスクは、アメフトにおける1つのプレーの中で、ボールを持っている選手(Rusher)が獲得したヤード数を予測するというものでした。

kernelでデータを確認すると、22選手(オフェンス、ディフェンス11人ずつ)のプレー開始時の位置や速度などプレー情報のほか、選手に紐づく年齢などの情報もデータフレームとして与えらていることがわかります。

選手の位置や動きをplotして可視化してしたものがsolutionでも共有されているので見てみると、イメージがつかめると思います。

こうしてデータを見ているとRusherの近くにディフェンスがいるかどうかなど、位置情報がとても重要そうだということがわかります。

1st place solutionのユニークさ

この問題を解くにあたり、解法はいくつか考えられると思います。

Rusherとの距離などを踏まえたhand craftな特徴を作って、BoostingやNN
GCNやTransformerを使って、上記の位置情報などの表現をNNで獲得

今回のコンペの上位は基本的に2の戦略を取っていたように思えますが、1st solutionはほとんど特徴を作らずにCNNで学習していて、他チームを大きく離すようなスコアを出していました。

f:id:to78314910:20191213005200p:plain — モデル：https://www.kaggle.com/c/nfl-big-data-bowl-2020/discussion/119400#latest-689317

こちらの画像もdiscussionからの引用で、1st solutionの全体像を表しています。流れは下記のようになります。

[バッチサイズ, 特徴, オフェンス, ディフェンス]のshapeを持ったバッチを作る
1つ目のCNNでオフェンス特徴を獲得する
2つ目のCNNでディフェンス特徴を抽出する
最後に獲得した特徴をTargetに変換する

特徴は5種類でディフェンスの速度、ディフェンスとRusherの位置・速度の距離、オフェンスとディフェンスの位置・速度の距離となっています。

再現チャレンジ

ここからはsolutionを、自分が再現実装したコードと照らし合わせてみていきます。

前処理

Data preprocessingのところで以下のような説明がありました。

targetを元の-99~99の範囲から、-30~50の範囲に修正。 SをDis*10で置き換える。(元のSの値に問題がある) フィールドの左向きのプレーと右向きのプレーがあるので、向きを統一する。

向きの統一はCPMPさんが公開してくれているkernelからreorient関数を持ってきて、その他は以下のように実装しました。

y_mae = train[TARGET_COLUMNS][::22].values
y_mae = np.where(y_mae < -30, -30, y_mae)
y_mae = np.where(y_mae > 50, 50, y_mae)

df["S"] = df["Dis"] * 10

特徴

特徴は先述の通り5種類あり以下のように実装しました。

ディフェンスの速度 -> X, Y方向の速度をそのまま使う
ディフェンスとRusherの位置、速度の距離 -> 位置は差分の絶対値、速度は差分をそのまま使う
オフェンスとディフェンスの位置、速度の距離 -> 上と同じ

特徴については絶対値にするか否かを色々実験して、local scoreを見る限り、これが最適かなという結果でした。X, Yを合わせてユークリッド距離を求めて特徴を得るという方法もあると思いますが、The Zooのメンバーに直接聞いたところ10特徴だと言っていました。

また、ここの実装は中国のKaggler、903124さんにも手伝ってもらいました。

dist_def_off_x = abs(def_x.reshape(-1,1)-off_x.reshape(1,-1))
dist_def_off_sx = def_sx.reshape(-1,1)-off_sx.reshape(1,-1)
dist_def_off_y = abs(def_y.reshape(-1,1)-off_y.reshape(1,-1))
dist_def_off_sy = def_sy.reshape(-1,1)-off_sy.reshape(1,-1)
dist_def_rush_x = abs(def_x.reshape(-1,1)-np.repeat(xysdir_rush[0],10).reshape(1,-1))
dist_def_rush_y = abs(def_y.reshape(-1,1)-np.repeat(xysdir_rush[1],10).reshape(1,-1))
dist_def_rush_sx = def_sx.reshape(-1,1)-np.repeat(xysdir_rush[2],10).reshape(1,-1)
dist_def_rush_sy = def_sy.reshape(-1,1)-np.repeat(xysdir_rush[3],10).reshape(1,-1)
def_sx = np.repeat(def_sx,10).reshape(11,-1)
def_sy = np.repeat(def_sy,10).reshape(11,-1)

スケーリング

スケーリングについてはあまり記載がなかったですが、チャンネルごとに特徴が異なるので、分布も異なるだろうということで、チャンネル毎のスケーリングを実装しました。

def scaling(feats, sctype="standard"):
    v1 = []
    v2 = []
    for i in range(feats.shape[1]):
        feats_ = feats[:, i, :]
        if sctype == "standard":
            mean_ = np.mean(feats_)
            std_ = np.std(feats_)
            feats[:, i, :] -= mean_
            feats[:, i, :] /= std_
            v1.append(mean_)
            v2.append(std_)
        elif sctype == "minmax":
            max_ = np.max(feats_)
            min_ = np.min(feats_)
            feats[:, i, :] = (feats_ - min_) / (max_ - min_)
            v1.append(max_)
            v2.append(min_)

    return feats, v1, v2

一応MinMaxScalingも試しましたが、StandardScalingの方がよかったです。

モデル

実装は以下のようにlayer構成を模倣しました。

最初のConv×3でチャンネル方向に特徴の交互作用のみを獲得するlayerを作る
次のPoolingでオフェンス方向に獲得した特徴をまとめる
次のConv×3で畳み込んだ特徴の交互作用のみを獲得するlayerを作る
最後のPoolingで獲得した特徴をディフェンス方向にまとめる
全結合層で獲得した特徴をTargetへ変換する

class CnnModel(nn.Module):

    def __init__(self, num_classes):

        super().__init__()

        self.conv1 = nn.Sequential(

            nn.Conv2d(10, 128, kernel_size=1, stride=1, bias=False),

            nn.ReLU(inplace=True),

            nn.Conv2d(128, 160, kernel_size=1, stride=1, bias=False),

            nn.ReLU(inplace=True),

            nn.Conv2d(160, 128, kernel_size=1, stride=1, bias=False),

            nn.ReLU(inplace=True)

        )

        self.pool1 = nn.AdaptiveAvgPool2d((1, 11))

 

        self.conv2 = nn.Sequential(

            nn.BatchNorm2d(128),

            nn.Conv2d(128, 160, kernel_size=(1, 1), stride=1, bias=False),

            nn.ReLU(inplace=True),

            nn.BatchNorm2d(160),

            nn.Conv2d(160, 96, kernel_size=(1, 1), stride=1, bias=False),

            nn.ReLU(inplace=True),

            nn.BatchNorm2d(96),

            nn.Conv2d(96, 96, kernel_size=(1, 1), stride=1, bias=False),

            nn.ReLU(inplace=True),

            nn.BatchNorm2d(96),

        )

        self.pool2 = nn.AdaptiveAvgPool2d((1, 1))

 

        self.last_linear = nn.Sequential(

            Flatten(),

            nn.Linear(96, 256),

            nn.LayerNorm(256),

            nn.ReLU(inplace=True),

            nn.Dropout(0.3),

            nn.Linear(256, num_classes)

        )

 

    def forward(self, x):

        x = self.conv1(x)

        x = self.pool1(x)

        x = self.conv2(x)

        x = self.pool2(x)

        x = self.last_linear(x)

 

        return x

1stの実装ではPoolingをGAP×0.7+GMP×0.3にしていましたが、そこをGAPのみに変更してあります。ただ、大きな違いはないと思います。

バリデーション

バリデーションについては、テストデータに当たる2019年のデータが、学習データのうち2017年より2018年のデータの傾向に近かったため、以下のような戦略を使用していると言っていました。

2018年データのみでlocalのscoreは計算
GameIDで5分割のGroupKFold
2017年データは学習のみに使用

そこで、下記のように実装しました。

x_2017, y_crps_2017, y_mae_2017 = x[season==2017], y_crps[season==2017], y_mae[season==2017]

x_usage, y_crps_usage, y_mae_usage = x[season!=2017], y_crps[season!=2017], y_mae[season!=2017]

folds = GroupKFold(n_splits=N_SPLITS).split(y_mae_usage, y_mae_usage, groups=game_id[season!=2017])

 

for n_fold, (train_idx, val_idx) in enumerate(folds):

    x_train, y_train = x_usage[train_idx], y_crps_usage[train_idx]

    x_val, y_val = x_usage[val_idx], y_crps_usage[val_idx]


    # add 2017 data

    x_train = np.concatenate([x_train, x_2017], axis=0)

    y_train = np.concatenate([y_train, y_crps_2017], axis=0)

Augmentation

アメフトのフィールド構造についてはコンペのデータページを見るとよくわかります。

特徴的なのは、アメフトのゲームにおいて、Y方向にデータが反転していても全く同じ結果になるだろうという点です。なので、Y軸方向のFlipを行うことができると、ソリューションでも述べられています。このAugmentationは学習時にもテスト時にも利用できます。

この部分についてはYに合わせてDirection(動く方向)も変える必要があるということで、そこがまだ反映できませんでした。

その他実装したこと

OptimizerとScheduler、Lossは下記のように実装しましたが、lossとSchedulerは自分の実装ではworkしていません。 OneCycleLRはこちらのgithubから転用しています。

class CRPSLoss(nn.Module):

    def __init__(self, n_class=199):

        super().__init__()

        self.n_class = n_class

 

    def forward(self, y_pred, y_true):

        y_true = torch.clamp(torch.cumsum(y_true, 1), 0, 1)

        y_pred = torch.clamp(torch.cumsum(y_pred, 1), 0, 1)

        crps = torch.mean((y_true - y_pred) ** 2)

        return crps

 

criterion = CRPSLoss()

optimizer = torch.optim.Adam(model.parameters(), lr=1e-3)

scheduler = OneCycleLR(optimizer, num_steps=num_steps, lr_range=(0.0005, 0.001))

結果

1stは最終的にlocalで0.01215だったと言っていましたが、自分の手元では0.01315程度でした。。。

Augumentationと4倍のseed averagingが実装できていないとは言え低い、、、

コードの公開が待たれますが、ルールとかもあるので、できるかわからないとのことでした。再現実装のチームメイト募集中です！