xonstat/glicko.py

   1 import logging
   2 import math
   3 import sys
   4
   5 log = logging.getLogger(__name__)
   6
   7 # DEBUG
   8 # log.addHandler(logging.StreamHandler())
   9 # log.setLevel(logging.DEBUG)
  10
  11 # the default initial rating value
  12 MU = 1500
  13
  14 # the default ratings deviation value
  15 PHI = 350
  16
  17 # the default volatility value
  18 SIGMA = 0.06
  19
  20 # the default system volatility constant
  21 #TAU = 0.3
  22 TAU = 0.5
  23
  24 # the ratio to convert from/to glicko2
  25 GLICKO2_SCALE = 173.7178
  26
  27
  28 class PlayerGlicko(object):
  29     def __init__(self, mu=MU, phi=PHI, sigma=SIGMA):
  30         self.mu = mu
  31         self.phi = phi
  32         self.sigma = sigma
  33
  34     def to_glicko2(self):
  35         """ Convert a rating to the Glicko2 scale. """
  36         return PlayerGlicko(
  37             mu=(self.mu - MU) / GLICKO2_SCALE,
  38             phi=self.phi / GLICKO2_SCALE,
  39             sigma=self.sigma
  40         )
  41
  42     def from_glicko2(self):
  43         """ Convert a rating to the original Glicko scale. """
  44         return PlayerGlicko(
  45             mu=self.mu * GLICKO2_SCALE + MU,
  46             phi=self.phi * GLICKO2_SCALE,
  47             sigma=self.sigma
  48         )
  49
  50
  51 def calc_g(phi):
  52     return 1 / math.sqrt(1 + (3 * phi ** 2) / (math.pi ** 2))
  53
  54
  55 def calc_e(mu, mu_j, phi_j):
  56     return 1. / (1 + math.exp(-calc_g(phi_j) * (mu - mu_j)))
  57
  58
  59 def calc_v(gs, es):
  60     """ Estimated variance of the team or player's ratings based only on game outcomes. """
  61     total = 0.0
  62     for i in range(len(gs)):
  63         total += (gs[i] ** 2) * es[i] * (1-es[i])
  64
  65     return 1. / total
  66
  67
  68 def calc_delta(v, gs, es, results):
  69     """
  70     Compute the estimated improvement in rating by comparing the pre-period rating to the
  71     performance rating based only on game outcomes.
  72     """
  73     total = 0.0
  74     for i in range(len(gs)):
  75         total += gs[i] * (results[i] - es[i])
  76
  77     return v * total
  78
  79
  80 def calc_sigma_bar(sigma, delta, phi, v, tau=TAU):
  81     """ Compute the new volatility. """
  82     epsilon = 0.000001
  83     A = a = math.log(sigma**2)
  84
  85     # pre-compute some terms
  86     delta_sq = delta ** 2
  87     phi_sq = phi ** 2
  88
  89     def f(x):
  90         e_up_x = math.e ** x
  91         term_a = (e_up_x * (delta_sq - phi_sq - v - e_up_x)) / (2 * (phi_sq + v + e_up_x) ** 2)
  92         term_b = (x - a) / tau ** 2
  93         return term_a - term_b
  94
  95     if delta_sq > (phi_sq + v):
  96         B = math.log(delta_sq - phi_sq - v)
  97     else:
  98         k = 1
  99         while f(a - k * tau) < 0:
 100             k += 1
 101         B = a - k * tau
 102
 103     fa, fb = f(A), f(B)
 104     while abs(B - A) > epsilon:
 105         C = A + (A - B) * (fa / (fb - fa))
 106         fc = f(C)
 107
 108         if fc * fb < 0:
 109             A, fa = B, fb
 110         else:
 111             fa /= 2
 112
 113         B, fb = C, fc
 114
 115         log.debug("A={}, B={}, C={}, fA={}, fB={}, fC={}".format(A, B, C, fa, fb, fc))
 116
 117     return math.e ** (A / 2)
 118
 119
 120 def rate(player, opponents, results):
 121     """
 122     Calculate the ratings improvement for a given player, provided their opponents and
 123     corresponding results versus them.
 124     """
 125     p_g2 = player.to_glicko2()
 126
 127     gs = []
 128     es = []
 129     for i in range(len(opponents)):
 130         o_g2 = opponents[i].to_glicko2()
 131         gs.append(calc_g(o_g2.phi))
 132         es.append(calc_e(p_g2.mu, o_g2.mu, o_g2.phi))
 133
 134         # DEBUG
 135         # log.debug("j={} muj={} phij={} g={} e={} s={}"
 136                   # .format(i+1, o_g2.mu, o_g2.phi, gs[i], es[i], results[i]))
 137
 138     v = calc_v(gs, es)
 139     delta = calc_delta(v, gs, es, results)
 140     sigma_bar = calc_sigma_bar(p_g2.sigma, delta, p_g2.phi, v)
 141
 142     phi_tmp = math.sqrt(p_g2.phi ** 2 + sigma_bar ** 2)
 143     phi_bar = 1/math.sqrt((1/phi_tmp**2) + (1/v))
 144
 145     sum_terms = 0.0
 146     for i in range(len(opponents)):
 147         sum_terms += gs[i] * (results[i] - es[i])
 148
 149     mu_bar = p_g2.mu + phi_bar**2 * sum_terms
 150
 151     new_rating = PlayerGlicko(mu_bar, phi_bar, sigma_bar).from_glicko2()
 152
 153     # DEBUG
 154     # log.debug("v={}".format(v))
 155     # log.debug("delta={}".format(delta))
 156     # log.debug("sigma_temp={}".format(sigma_temp))
 157     # log.debug("sigma_bar={}".format(sigma_bar))
 158     # log.debug("phi_bar={}".format(phi_bar))
 159     # log.debug("mu_bar={}".format(mu_bar))
 160     # log.debug("new_rating: {} {} {}".format(new_rating.mu, new_rating.phi, new_rating.sigma))
 161
 162     return new_rating
 163
 164
 165 def main():
 166     pA = PlayerGlicko(mu=1500, phi=200)
 167     pB = PlayerGlicko(mu=1400, phi=30)
 168     pC = PlayerGlicko(mu=1550, phi=100)
 169     pD = PlayerGlicko(mu=1700, phi=300)
 170
 171     opponents = [pB, pC, pD]
 172     results = [1, 0, 0]
 173
 174     rate(pA, opponents, results)
 175
 176
 177 if __name__ == "__main__":
 178     sys.exit(main())