箱根駅伝で優勝する確率をエントリー時のタイムを用いて予測してみる

Stan python

はじめに

箱根駅伝は好きで割と毎年見ているのですが、 マラソン系の競技ってパフォーマンスの分散が対人競技と比較すると小さくなるため、実は事前の情報から結果を予測することが比較的容易なんじゃないか？という疑問が浮かんだので、各大学が優勝する確率を出してみることにしました。

方針としては、箱根駅伝に出場する各校のトップ10人の10000m平均タイムより箱根駅伝の予測タイムの分布を算出し、 算出した予測分布を用いて各校が優勝する確率を予測します。

データの収集・確認

各年の10000mのタイムと箱根駅伝の総合タイムについては以下のサイト参考にさせていただきました。

上記のサイトより2013-2018年の過去５年分の箱根駅伝エントリー時の各校上位10人の10000m平均タイムを取得し、 ２次元プロットした図が以下になります。

f:id:rmizutaa:20190107183326p:plain

横軸が10000mの平均タイム、縦軸が箱根駅伝の総合タイムです。単位は秒に直しています。 2つの変数には正の相関がありますが、それなりにばらつきもありそうです。

モデルの構築

今回求めたいのは予測値ではなく予測分布なので、機械学習ではなく統計モデリングの手法を用います。 具体的には、ある統計モデルを作成し、そのモデルのパラメータを推定するためにStanでMCMCを行います。

統計モデリングの詳しい話は割愛しますが、 StanとRでベイズ統計モデリング とかがわかりやすいみたいです（まだ読了してません。。。）

今回はモデルを以下のようにおきました。

y = normal(a*X[N], \sigma)

yを箱根駅伝の総合タイム、X[N]をNサンプル目の10000mタイムとし、 yが平均a*X[N]、分散\sigmaの正規分布であると仮定しています。 stanを用いてパラメータaと\sigmaの推定を行います。

stanのコードは以下のようになります。Rは不得手なのでpystanを使っています。

model = """
    data {
        int<lower=0> N;
        real X[N];
        real Y[N];
    }
    parameters {
        real a;
        real<lower=0> sigma;
    }
    model { 
        for (i in 1:N)
            Y[i] ~ normal(a*X[N], sigma);
    }
"""

stanを実行します。

stan_data = {'N': len(x), 'X':x, 'Y': y}
sm = pystan.StanModel(model_code=model)
fit = sm.sampling(data=stan_data, iter=10000, chains=3)

結果は下記のようになりました。 rhatはどの値も1.0なので収束していますね。

f:id:rmizutaa:20190107190812p:plain

勝率の算出

前章で導出したパラメータを使って各校の勝率を算出してみます。 2019年の10000m平均タイムと先ほど推定したパラメータを用いて各校の予測分布を作成し、 そこからランダムサンプリングを10000回行い、優勝する確率を予測します。

import numpy as np
from numpy.random import *

#試行回数
num=10000

#空箱
arr=np.zeros((num,len(time_10000m)))

#ランダムサンプリングにより各校の予測タイムを出力
for n in range(num):
    for order, ti in enumerate(time_10000m):
        arr[n,order]=normal(22.93*ti,623.67) #パラメータは先ほど求めた値
        
#各試行のトップの大学をカウント、割合に
prob=pd.DataFrame(pd.Series(arr.argmin(axis=1)).value_counts()/num)

result=df2018.join(prob)
result=result.sort_values(by=0,ascending=False)

#グラフ用
plt.figure(figsize=(12,8))
plt.bar(result["大学名"],result[0])
plt.xticks(rotation=30)
plt.show()

f:id:rmizutaa:20190108092548p:plain

今回の仮定の置き方だと優勝する確率は青山学院大が19%、東海大が12.5%、駒沢大、帝京大が9%、東洋大が7％という結果になりました。 (実際の順位は1位：東海大、2位：青山学院大、3位：東洋大、4位：駒沢大、5位：帝京大)

考察

個人的には10000m平均タイムが上位のチームの優勝確率はもう少し高く、下位のチームの優勝確率はもう少し低いイメージなので、実際より予測分布が広がってしまっていたかなという気がします。 原因としては、今回は上位10人の10000mの平均タイムのみから箱根駅伝の総合タイムを予測するという結構大雑把な仮定で行ったためだと考えられるので、天候や選手層の厚み、走者の学年分布等もモデルに追加することで、もう少し分散の小さい予測が期待できるかなと思います。

使用したコードは以下になります。
https://github.com/rmizuta3/hakone_ekiden