Main Code for "Features in Trueskill"

mcmc_ranking.py

import numpy as np

import theano.tensor as T
import pymc3 as pm

from utils import add_params_property

import ranking
from ranking import tennis_data, MPTrueSkill1V1NoDrawRanker

from sklearn.model_selection import cross_val_score

import scipy.stats


def logistic(x):
    return 1.0 / (1.0 + T.exp(-x))


def approximate_normal_cdf(z):
    # http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.429.6900&rep=rep1&type=pdf
    y = 0.07056*z**3 + 1.5976*z
    return logistic(y)


class WinRatio1V1Ranker:
    @add_params_property
    def __init__(
            self,
            n_players,
    ):
        self.params = vars(self.params)
        self.__dict__ = {**self.__dict__, **self.params}
        self.fitted = False

    def get_params(self, deep):
        return self.params

    def set_params(self, params):
        self.params = params

    def fit(self, match_data, outcomes):
        matches = match_data[:, :2]
        assert outcomes.shape == (len(matches), 3)
        n_wins = np.zeros(self.n_players)
        n_losses = np.zeros(self.n_players)

        for (p0, p1), outcome in zip(matches, outcomes):
            if outcome[0]:
                # p0 wins
                n_wins[p0] += 1
                n_losses[p1] += 1
            if outcome[2]:
                # p1 wins
                n_wins[p1] += 1
                n_losses[p0] += 1

        self.win_ratios = n_wins / (n_losses + 1)
        self.fitted = True

    def predict(self, match_data):
        matches = match_data[:, :2]
        p0_wins = self.win_ratios[matches[:, 0]
                                  ] >= self.win_ratios[matches[:, 1]]
        outcomes = np.zeros([len(match_data), 3], dtype=bool)
        outcomes[:, 0] = p0_wins
        outcomes[:, 2] = np.logical_not(p0_wins)
        return outcomes


class MCMCTrueSkill1V1Ranker:
    @add_params_property
    def __init__(
            self,
            n_players,
            prior_env_std,
            prior_skill_std,
            prior_draw_eps,
    ):
        self.params = vars(self.params)
        self.__dict__ = {**self.__dict__, **self.params}
        self.fitted = False

    def get_params(self, deep):
        return self.params

    def set_params(self, params):
        self.params = params

    def fit(self, match_data, outcomes_observed):
        matches = match_data[:, :2]
        print('Fitting')
        assert outcomes_observed.shape == (len(matches), 3)

        with pm.Model() as model:
            skills = pm.Normal(
                'skills', mu=0, sd=self.prior_skill_std, shape=self.n_players)
            skill_differences = skills[matches[:, 0]
                                       ] - skills[matches[:, 1]]

            sd = self.prior_env_std
            skill_differences = skill_differences
            draw_eps = self.prior_draw_eps
            p_ds_lt_minus_eps = approximate_normal_cdf(
                (-draw_eps - skill_differences) / sd)
            p_ds_lt_plus_eps = approximate_normal_cdf(
                (+draw_eps - skill_differences) / sd)

            p1_win_probs = p_ds_lt_minus_eps
            draw_probs = p_ds_lt_plus_eps - p_ds_lt_minus_eps
            p0_win_probs = 1 - p_ds_lt_plus_eps

            outcome_probs = T.stack(
                [p0_win_probs, draw_probs, p1_win_probs], axis=-1)

            outcomes = pm.Categorical(
                'outcomes',
                p=outcome_probs,
                observed=np.argmax(outcomes_observed, axis=-1)
            )

            self.trace = pm.sample(1000)

        self.fitted = True

    def predict(self, match_data):
        matches = match_data[:, :2]
        skills = self.trace.get_values('skills', burn=500, thin=2)
        skill_d = skills[:, matches[:, 0]] - skills[:, matches[:, 1]]
        p0_wins_sim = skill_d > 0
        p0_wins_pred = np.mean(p0_wins_sim, axis=0) > 0.5
        preds = np.zeros([len(matches), 3], dtype=bool)
        preds[:, 0] = p0_wins_pred
        preds[:, 2] = np.logical_not(p0_wins_pred)
        return preds

    def predict(self, match_data):
        matches = match_data[:, :2]
        p0_wins = self.q_skills[matches[:, 0]] > self.q_skills[matches[:, 1]]
        outcomes = np.zeros([len(match_data), 3], dtype=bool)
        outcomes[:, 0] = p0_wins
        outcomes[:, 2] = np.logical_not(p0_wins)
        return outcomes


class MAPTrueSkill1V1RankerWithP0Advantage:
    @add_params_property
    def __init__(
            self,
            n_players,
            prior_env_std,
            prior_skill_std,
            prior_draw_eps,
            prior_p0_adv_std,
            advantage_prior,
    ):
        self.params = vars(self.params)
        self.__dict__ = {**self.__dict__, **self.params}
        self.fitted = False

    def get_params(self, deep):
        return self.params

    def set_params(self, params):
        self.params = params

    def fit(self, match_data, outcomes_observed):
        matches = match_data[:, :2]
        print('Fitting')
        assert outcomes_observed.shape == (len(matches), 3)
        if self.advantage_prior is None:
            p0_win_percentage = np.mean(outcomes_observed[:, 0])
            prior_win_p0_adv = \
                scipy.stats.norm.ppf(p0_win_percentage)*self.prior_env_std
        else:
            prior_win_p0_adv = self.advantage_prior

        print('p0 adv prior', prior_win_p0_adv)
        with pm.Model() as model:
            skills = pm.Normal(
                'skills', mu=0, sd=self.prior_skill_std, shape=self.n_players)
            skill_differences = skills[matches[:, 0]
                                       ] - skills[matches[:, 1]]

            p0_adv = pm.Normal('p0_adv', mu=prior_win_p0_adv, sd=self.prior_p0_adv_std)
            sd = self.prior_env_std
            skill_differences = skill_differences + p0_adv
            draw_eps = self.prior_draw_eps
            p_ds_lt_minus_eps = approximate_normal_cdf(
                (-draw_eps - skill_differences) / sd)
            p_ds_lt_plus_eps = approximate_normal_cdf(
                (+draw_eps - skill_differences) / sd)

            p1_win_probs = p_ds_lt_minus_eps
            draw_probs = p_ds_lt_plus_eps - p_ds_lt_minus_eps
            p0_win_probs = 1 - p_ds_lt_plus_eps

            outcome_probs = T.stack(
                [p0_win_probs, draw_probs, p1_win_probs], axis=-1)

            outcomes = pm.Categorical(
                'outcomes',
                p=outcome_probs,
                observed=np.argmax(outcomes_observed, axis=-1)
            )

            map_estimate = pm.find_MAP(maxeval=10000)
            self.q_skills = map_estimate['skills']
            self.p0_adv = map_estimate['p0_adv']

        self.fitted = True

    def predict(self, match_data):
        matches = match_data[:, :2]
        p0_wins = self.q_skills[matches[:, 0]] + \
            self.p0_adv > self.q_skills[matches[:, 1]]
        outcomes = np.zeros([len(match_data), 3], dtype=bool)
        outcomes[:, 0] = p0_wins
        outcomes[:, 2] = np.logical_not(p0_wins)
        return outcomes


class MAPTrueSkill1V1RankerWithAffinity:
    @add_params_property
    def __init__(
            self,
            n_players,
            prior_env_std,
            prior_skill_std,
            prior_draw_eps,
            feature_weight,
    ):
        self.params = vars(self.params)
        self.__dict__ = {**self.__dict__, **self.params}
        self.fitted = False

    def get_params(self, deep):
        return self.params

    def set_params(self, params):
        self.params = params

    def fit(self, match_data, outcomes_observed):
        print('Fitting')
        matches = match_data[:, :2]
        assert outcomes_observed.shape == (len(matches), 3)

        match_types = match_data[:, 2]
        ohe_match_types = np.zeros([len(match_types), 3])
        ohe_match_types[np.arange(len(match_types)), match_types] = 1.0

        with pm.Model() as model:
            skills = pm.Normal(
                'skills', mu=0, sd=self.prior_skill_std, shape=self.n_players)
            skill_differences = \
                skills[matches[:, 0]] - skills[matches[:, 1]]

            affinities = pm.Normal('affinities', mu=0,
                                   sd=self.prior_skill_std, shape=[self.n_players, 3])

            ms_affinities = affinities - \
                T.mean(affinities, axis=-1, keepdims=True)

            affin_differences = \
                ms_affinities[matches[:, 0]] - ms_affinities[matches[:, 1]]

            relevant_affin_differences = \
                T.sum((affin_differences * ohe_match_types), axis=-1)

            sd = self.prior_env_std
            a = self.feature_weight
            skill_differences = \
                (1-a)*skill_differences + a*relevant_affin_differences
            draw_eps = self.prior_draw_eps
            p_ds_lt_minus_eps = approximate_normal_cdf(
                (-draw_eps - skill_differences) / sd)
            p_ds_lt_plus_eps = approximate_normal_cdf(
                (+draw_eps - skill_differences) / sd)

            p1_win_probs = p_ds_lt_minus_eps
            draw_probs = p_ds_lt_plus_eps - p_ds_lt_minus_eps
            p0_win_probs = 1 - p_ds_lt_plus_eps

            outcome_probs = T.stack(
                [p0_win_probs, draw_probs, p1_win_probs], axis=-1)

            outcomes = pm.Categorical(
                'outcomes',
                p=outcome_probs,
                observed=np.argmax(outcomes_observed, axis=-1)
            )

            self.trace = pm.sample(1000)

        self.fitted = True

    def predict(self, match_data):
        matches = match_data[:, :2]

        match_types = match_data[:, 2]
        ohe_match_types = np.zeros([len(match_types), 3])
        ohe_match_types[np.arange(len(match_types)), match_types] = 1.0

        skills = self.trace.get_values('skills', burn=500, thin=2)
        affinities = self.trace.get_values('affinities', burn=500, thin=2)
        affinities -= np.mean(affinities, axis=-1, keepdims=True)
        all_aff_d = affinities[:, matches[:, 0]] - affinities[:, matches[:, 1]]
        aff_d = np.sum(all_aff_d * ohe_match_types, axis=-1)
        skill_d = skills[:, matches[:, 0]] - skills[:, matches[:, 1]]
        p0_wins_sim = aff_d + skill_d > 0
        p0_wins_pred = np.mean(p0_wins_sim, axis=0) > 0.5
        preds = np.zeros([len(matches), 3], dtype=bool)
        preds[:, 0] = p0_wins_pred
        preds[:, 2] = np.logical_not(p0_wins_pred)
        return preds


def main():

    datasets = [
        ['ATP Tour', ranking.tennis_data()],
        #  ['NFL', ranking.nfl_data()],
        #  ['Online Chess', ranking.chess_data()],
    ]
    for _ in range(4):
        synth = [
            #  ['Synthetic', ranking.synthetic_data(
            #      n_players=100,
            #      n_matches=500,
            #      skill_std=1.0,
            #      env_noise_std=1.0,
            #      draw_eps=0.1,
            #  )],
            #  ['Synthetic Adv', ranking.synthetic_data_with_p0_advantage(
            #      n_players=100,
            #      n_matches=500,
            #      skill_std=1.0,
            #      env_noise_std=0.5,
            #      draw_eps=0.1,
            #      percent_equal_skill_p0_wins=0.85,
            #  )],
            ['Synthetic Aff',
             ranking.synthetic_data_with_match_types(
                 skill_std=1.0,
                 env_noise_std=0.5,
                 draw_eps=0.0,
                 n_players=35,
                 n_matches=600,
             )],
        ]
        datasets += synth
    for dataset_name, dataset in datasets:
        match_data, outcomes, n_players = dataset
        not_draws = np.logical_not(outcomes[:, 1])
        match_data = match_data[not_draws]
        outcomes = outcomes[not_draws]
        models = [
            ['Win Ratio', WinRatio1V1Ranker(
                n_players=n_players,
            )],
            ['Trueskill - Message Passing', ranking.MPTrueSkill1V1NoDrawRanker(
                n_players=n_players,
            )],
            #  ['Trueskill - MAP', MAPTrueSkill1V1Ranker(
            #      n_players=n_players,
            #      prior_env_std=np.sqrt(2),
            #      prior_skill_std=1.0,
            #      prior_draw_eps=0.001,
            #  )],
            #  ['Trueskill with P0 Advantage - MAP', MAPTrueSkill1V1RankerWithP0Advantage(
            #      n_players=n_players,
            #      prior_env_std=np.sqrt(2),
            #      prior_skill_std=1.0,
            #      prior_draw_eps=0.001,
            #      prior_p0_adv_std=1.0,
            #      advantage_prior=(0.0 if 'Chess' in dataset_name else None),
            #  )],
            ['Trueskill with Match Type Affinity - MAP',
                MAPTrueSkill1V1RankerWithAffinity(
                    n_players=n_players,
                    prior_env_std=np.sqrt(2),
                    prior_skill_std=1.0,
                    prior_draw_eps=0.001,
                    feature_weight=0.2,
                )
             ]
        ]
        for model_name, model in models:
            if 'Affinity' in model_name and match_data.shape[1] != 3:
                continue
            cvs = cross_val_score(model, match_data, outcomes,
                                  scoring='accuracy', cv=10)
            results = {
                'model': model_name,
                'dataset': dataset_name,
                'cvs': list(cvs),
                'mean': np.mean(cvs),
                '+2std': np.mean(cvs)+2*np.std(cvs),
                '-2std': np.mean(cvs)-2*np.std(cvs),
            }
            print(results)


if __name__ == '__main__':
    main()