<p><div><div><div><div data-show-faces="false" data-height="21" data-width="100" data-layout="button_count" data-send="true" data-href="https://www.r-bloggers.com/elo-and-elobeta-models-in-snooker/"></div></div><div></div></div><div></div></div>
<div>
(This article was first published on <strong><a href="http://www.questionflow.org/2018/07/03/elo-and-elobeta-models-in-snooker/"> QuestionFlow </a></strong>, and kindly contributed to <a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers)</a> </div>
<p> Research about adequacy of Elo based models applied to snooker match results. Contains a novel approach (EloBeta) targeted for sport results with variable “best of N” format.<br>
</p>
<div>
<h1>Prologue</h1>
<p>For many years I’ve been following snooker as a sport. It has it all: hypnotic beauty of smart play, elegance of cue strikes and psychological tension of competition. The one thing I don’t like that much is its <em>ranking system</em>.</p>
<p>Generally speaking, <strong>current snooker ranking is based on player accomplishments in tournaments (events)</strong> with different “weights” for different tournaments. Long time ago, it just used World Championships. Then, after more events had emerged, there was a table of points player could earn for winning at certain stage of tournament. Now it has the form of “rolling” sum of prize money player won during (roughly) past two calendar years.</p>
<p>This system has two main advantages: it is <strong>simple</strong> (win more money -> rise in rankings) and <strong>predictabile</strong> (want to get certain ranking -> win certain amount of money, other things being equal). The problem is that <strong>this type of rankings doesn’t account the strength (skill, form) of player’s opponents</strong>. The usual counter-argument for this is that if player reached high stage of tournament then he/she is “strong” at this moment of time (“weak players don’t win tournaments”). Well, it does sound quite convincing. However, in snooker, as in any sport, the role of chance should be taken into account: if player is “weaker” it doesn’t mean that he can’t ever win in a match with “stronger” player. It means that this happens less often then the other way around. Here where <a rel="nofollow" href="https://en.wikipedia.org/wiki/Elo_rating_system" target="_blank">Elo model</a> comes into play.</p>
<p>The idea behind Elo model is that each player is associated with numerical rating. The assumption is that a result of a game between two players can be predicted based on difference of their ratings: more value indicates more probability of “stronger” (with higher rating) player to win. <strong>Elo ranking is based on current player “strength”</strong> derived by wins against other players. This eliminates main disadvantage of current official ranking system. It is also capable of updating player rating during tournament to numerically react to player’s strong tournament performance.</p>
<p>Having some practical experience with Elo ratings, I think it can do well in snooker too. However, there is one obstacle: <em>it is devised for competitions with uniform type of games</em>. Yes, there are some variations to account for <a rel="nofollow" href="https://medium.com/@mattbarger/soccer-elo-the-rebuild-df6b58bd8b94" target="_blank">home field advantage</a> in football or <a rel="nofollow" href="https://arxiv.org/pdf/1012.4571.pdf" target="_blank">first move advantage</a> in chess (both by adding fixed amount of rating points to “less advantageous” player). In snooker, however, matches are played in the “best of <span>\(N\)</span>” format: the first one to win <span>\(n = \frac{N + 1}{2}\)</span> frames wins a match. We will also call this format “<span>\(n\)</span> to win”.</p>
<p>Intuitively, winning a “10 to win” match (final of major tournament) should be harder for “weaker” player then “4 to win” match (first rounds of current Home Nations tournaments). This is taken into account by my proposed <strong>EloBeta model</strong>.</p>
<p>To celebrate actual start of 2018/19 snooker season I decided to write this post in which <strong>I explore adequacy of both Elo and EloBeta models on snooker match results</strong>. Note that the goal is not to build models for forecasting and gambling purposes but for assessing players “strength” and creating “fair” ranking.</p>
<p>The idea of using Elo rating in snooker is not new at all. There are works on this topic, for example:</p>
<ul>
<li><a rel="nofollow" href="http://www.snookeranalyst.com/current-ratings/rating-summary/" target="_blank">Snooker Analyst</a> provides “Elo inspired” (more like <a rel="nofollow" href="https://en.wikipedia.org/wiki/Bradley-Terry_model" target="_blank">Bradley–Terry model</a>) rating system based on the idea of updating rating based on difference between “actual frames won” and “expected frames won”. This approach is a little bit questionable. Surely, more frame difference should indicate more difference in strength, however, achieving that is not player’s goal. In snooker aim is “just” to win match, i.e. get certain amount of frame wins before the opponent.</li>
<li><a rel="nofollow" href="https://www.snookerisland.com/forum/viewtopic.php?f=59&t=5585" target="_blank">This forum discussion</a> with implementation of basic Elo model.</li>
<li>Possibly, there are other works that I’ve missed. I will highly appreciate any information on this topic.</li>
</ul>
</div>
<div>
<h1>Overview</h1>
<p>This post is intended for both R users interested in Elo ratings and snooker analysts and fans. All experiments are written with intention to be reproducible. All code is hidden under spoilers (text appears after clicking on its summary, usually beginning with “Code for …”). It has commentaries and uses <a rel="nofollow" href="https://www.tidyverse.org/" target="_blank">tidyverse</a> packages, so it might be an interesting read for R users and programmers.</p>
<p>This post is organized as follows:</p>
<ul>
<li><strong>Models</strong> describes <em>Elo</em> and <em>EloBeta</em> models with their R implementations.</li>
<li><strong>Experiment</strong> describes details and intentions of computations: which <em>data</em> and <em>methodology</em> were used (and why) and what are the <em>results</em>.</li>
<li><strong>Exploration of EloBeta ratings</strong> has application results of EloBeta model to actual snooker data. This section is written more for Snooker fans than for R enthusiasts.</li>
</ul>
<p>We will need the following setup:</p>
<details>

<div>
<pre><code># Data wrangling packages
suppressPackageStartupMessages(library(dplyr))
library(tidyr)
library(purrr)
# Visualization package
library(ggplot2)
# Package for ratings
suppressPackageStartupMessages(library(comperank))

theme_set(theme_bw())

# Shouldn't be needed. Added just in case.
set.seed(20180703)</code></pre>
</div>
</details>
</div>
<div>
<h1>Models</h1>
<p>Both models are based on the following assumptions:</p>
<ol>
<li>There is a fixed set of players which should be ranked from “strongest” (first place) to “weakest” (last place).</li>
<li>Ranking is done by associating player <span>\(i\)</span> with numerical rating <span>\(r_i\)</span>: a number indicating the “strength” of player (more value -> “stronger”).</li>
<li>The more difference in player ratings before the match the less favorable is “weaker” player to win it.</li>
<li>Ratings are updated after every match based on its result and the ratings before it.</li>
<li>Winning against “stronger” opponent should lead to bigger increase in rating than winning against “weaker” opponent. The opposite should be true for losing.</li>
</ol>
<div>
<h2>Elo</h2>
<details>

<div>
<pre><code>#' @details This function is vectorized by all its arguments. Also usage of
#' `...` is crucial to allow supplying unrelated arguments in the future.
#' 
#' @return A probability of player 1 (rating `rating1`) wins in a match with
#'   player 2 (rating `rating2`). Here difference in ratings directly affects
#'   the outcome.
elo_win_prob &lt;- function(rating1, rating2, ksi = 400, ...) {
  norm_rating_diff &lt;- (rating2 - rating1) / ksi
  
  1 / (1 + 10^norm_rating_diff)
}

#' @return A rating function for Elo model that can be supplied to
#'   `comperank::add_iterative_ratings()`.
elo_fun_gen &lt;- function(K, ksi = 400) {
  function(rating1, score1, rating2, score2) {
    comperank::elo(rating1, score1, rating2, score2, K = K, ksi = ksi)[1, ]
  }
}</code></pre>
</div>
</details>
<p>Elo model updates ratings by the following steps:</p>
<ul>
<li>
<p><strong>Compute probability</strong> (<em>before</em> the match) <strong>of certain player to win the match</strong>. Probability of one player (we will call him/her “first”) with “identifier” <span>\(i\)</span> and rating <span>\(r_i\)</span> winning against the other player (“second”) with “identifier” <span>\(j\)</span> and rating <span>\(r_j\)</span> is equal to</p>
<p><span>\[Pr(r_i , r_j) = \frac{1}{1 + 10^{(r_j – r_i)/400}}\]</span></p>
<p>This way of computing probability is aligned with third model assumption.</p>
<p>Difference normalization by 400 is a mathematical way to say which difference is considered “big”. This can be replaced by a model parameter <span>\(\xi\)</span>, however this only affects the spread of future ratings and is often an overkill. Number 400 is fairly standard in chess world.</p>
<p>In general approach probability is equal to <span>\(L(r_j – r_i)\)</span> where <span>\(L(x)\)</span> is some strictly increasing function with values from 0 to 1. We will use logistic curve to compute winning probability. More thorough study can be found in <a rel="nofollow" href="https://www.stat.berkeley.edu/~aldous/Papers/me-Elo-SS.pdf" target="_blank">this article</a>.</p>
</li>
<li><strong>Obtain match result</strong> <span>\(S\)</span>. In basic model it is 1 if first player wins (second player loses), 0.5 in case of draw and 0 if second player wins.</li>
<li><strong>Update ratings</strong>:
<ul>
<li><span>\(\delta = K \cdot (S – Pr(r_i , r_j))\)</span>. This is the value by which ratings will change. It introduces the “K factor” (only model parameter). Less <span>\(K\)</span> (ratings being equal) means less change in ratings – model is more conservative, i.e. more wins is needed to “prove” increase in “strength”. On the other hand, more <span>\(K\)</span> means more “trust” to the current results than current ratings. <em>Choosing “optimal” <span>\(K\)</span> is a way to produce “good” ranking system</em>.</li>
<li><span>\(r_i^{(new)} = r_i + \delta\)</span>, <span>\(r_j^{(new)} = r_j – \delta\)</span>.</li>
</ul>
</li>
</ul>
<p><strong>Notes</strong>:</p>
<ul>
<li>As one can see from rating update formulas, the sum of ratings for all ranked players doesn’t change over time: rating increase of one rating can be only done by taking this amount from another player.</li>
<li>Players without any matches played are associated with initial rating 0. The usual value is 1500, however I don’t see any other reason except psychological for this. With previous note, using 0 means that sum of all ratings will always be 0, which is kind of beautiful.</li>
<li>It is needed some matches to be played in order for rating to represent player’s “strength”. This introduces a problem: newly added players start with rating 0 which is almost surely not the lowest among current players. In other words, newcomers are considered to be “stronger” than some other players. This should be dealt with by external procedures of rating updates when introducing new player: maybe he/she should start with some low rating while compensating overall sum decrease by uniform increasing of other players’ ratings.</li>
<li>Why this procedure makes sense? In case of equal ratings <span>\(\delta\)</span> always equals <span>\(0.5 \cdot K\)</span>. Let’s say, for example, that <span>\(r_i = 0\)</span> and <span>\(r_j = 400\)</span>. It means that probability of first player winning is <span>\(\frac{1}{1 + 10} \approx 0.0909\)</span>, i.e. he/she will win 1 out of 11 matches.
<ul>
<li>In case of win, he/she will be awarded with approximately <span>\(0.909 \cdot K\)</span> increase, which is more than in case of equal ratings.</li>
<li>If he/she is defeated, then rating is decreased by approximately <span>\(0.0909 \cdot K\)</span>, which is less than in case of equal ratings.</li>
</ul>
<p>This shows that Elo model is aligned with fifth model assumption: winning against “stronger” opponent leads to bigger increase in rating than winning against “weaker” opponent and vice versa.

<p>Of course, Elo model has its (fairly high-level) <a rel="nofollow" href="https://en.wikipedia.org/wiki/Elo_rating_system#Practical_issues" target="_blank">practical issues</a>. However, the most important for our research is that “it thinks” that all matches are played in uniform conditions. It means, that match length isn’t taken into account: winning in “4 to win” match is rewarded the same as winning in “10 to win” match. That is where EloBeta comes into play.</p>

<div>
<h2>EloBeta</h2>
<details>

<div>
<pre><code>#' @details This function is vectorized by all its arguments.
#' 
#' @return A probability of player 1 (rating `rating1`) wins in a match with
#'   player 2 (rating `rating2`). Match is played until either player wins
#'   `frames_to_win` frames. Here difference in ratings directly affects
#'   the probability of winning one frame.
elobeta_win_prob &lt;- function(rating1, rating2, frames_to_win, ksi = 400, ...) {
  prob_frame &lt;- elo_win_prob(rating1 = rating1, rating2 = rating2, ksi = ksi)
  
  # Probability that first player wins `frames_to_win` frames sooner than second
    # player based on probability of first player to win one frame `prob_frame`.
    # Frames are treated as independent games.
  pbeta(prob_frame, frames_to_win, frames_to_win)
}

#' @return Match result in terms of player 1 win: 1 if he/she wins, 0.5 in case
#'   of a draw, and 0 if he/she loses.
get_match_result &lt;- function(score1, score2) {
  # There are no ties in snooker but this handles general case
  near_score &lt;- dplyr::near(score1, score2)
  
  dplyr::if_else(near_score, 0.5, as.numeric(score1 > score2))
}

#' @return A rating function for EloBeta model that can be supplied to
#'   `comperank::add_iterative_ratings()`.
elobeta_fun_gen &lt;- function(K, ksi = 400) {
  function(rating1, score1, rating2, score2) {
    prob_win &lt;- elobeta_win_prob(
      rating1 = rating1, rating2 = rating2,
      frames_to_win = pmax(score1, score2), ksi = ksi
    )
    
    match_result &lt;- get_match_result(score1, score2)
    delta &lt;- K * (match_result - prob_win)
    
    c(rating1 + delta, rating2 - delta)
  }
}</code></pre>
</div>
</details>
<p>In Elo model difference in ratings directly affects the outcome probability of winning the whole match. The main idea behind EloBeta model is to <strong>make difference in ratings directly affect the outcome of one frame</strong> and <strong>actually compute the probability of a player winning <span>\(n\)</span> frames before his opponent</strong>.</p>
<p>The question is: how to compute this probability? Turns out, this is one of the oldest task in the history of probability theory and has its own name – <a rel="nofollow" href="https://en.wikipedia.org/wiki/Problem_of_points" target="_blank">Problem of points</a>. Very nice description can be found in <a rel="nofollow" href="https://probabilityandstats.wordpress.com/2016/11/06/the-problem-of-points/" target="_blank">this post</a>. Using its notation, the outcome probability equals to:</p>
<p><span>\[<br>
P(n, n) = \sum\limits_{j = n}^{2n-1}{{{2n-1} \choose {j}} p^j (1-p)^{2n-1-j}}<br>
\]</span></p>
<p>Here <span>\(P(n, n)\)</span> is a probability of first player to win in “<span>\(n\)</span> to win” match; <span>\(p\)</span> is his/her probability to win one frame (<span>\(1-p\)</span> – opponents probability). With this approach one assumes that <strong>result of frames inside a match is independent one from another</strong>. This is arguable, but necessary assumption (in terms of this model).</p>
<p>Is there a way to compute this faster? It turns out, that the answer is yes. After hours of formula wrangling, practical experiments and internet search I found the following <a rel="nofollow" href="https://dlmf.nist.gov/8.17#E5" target="_blank">property</a> of <a rel="nofollow" href="https://en.wikipedia.org/wiki/Beta_function#Incomplete_beta_function" target="_blank">regularized incomplete beta function</a> <span>\(I_x(a, b)\)</span>. By substituting <span>\(m = k,~ n = 2k – 1\)</span> into that property and changing <span>\(k\)</span> into <span>\(n\)</span> we obtain that <span>\(P(n, n) = I_p(n, n)\)</span>.</p>
<p>This is also very good news for R users, because <span>\(I_p(n, n)\)</span> can be computed as simply as <code>pbeta(p, n, n)</code>. <strong>Note</strong> that the general case probability of winning <span>\(n\)</span> frames before opponent wins <span>\(m\)</span> can also be computed as <span>\(I_p(n, m)\)</span> and <code>pbeta(p, n, m)</code> respectively. This introduces the reach field of <em>updating winning probabilities during the mach</em>.</p>
<p>So the procedure of updating ratings within EloBeta model is as follows (given ratings <span>\(r_i, r_j\)</span>, number of frames <span>\(n\)</span> to win in match and match outcome <span>\(S\)</span>, as in Elo model):</p>
<ul>
<li><strong>Compute probability of first player to win the frame</strong>: <span>\(p = Pr(r_i , r_j) = \frac{1}{1 + 10^{(r_j – r_i)/400}}\)</span>.</li>
<li><strong>Compute probability of this player to win the match</strong>: <span>\(Pr^{Beta}(r_i, r_j) = I_p(n, n)\)</span>. For example, if <span>\(p\)</span> equals 0.4 then probability to win “4 to win” match drops to 0.29 and with “18 to win” to 0.11.</li>
<li><strong>Update ratings</strong>:
<ul>
<li><span>\(\delta = K \cdot (S – Pr^{Beta}(r_i , r_j))\)</span>.</li>
<li><span>\(r_i^{(new)} = r_i + \delta\)</span>, <span>\(r_j^{(new)} = r_j – \delta\)</span>.</li>
</ul>
</li>
</ul>
<p><strong>Note</strong> that, as difference in ratings affects probability of winning one frame (not the whole match), we should expect lower optimal <span>\(K\)</span>: part of <span>\(\delta\)</span>’s value comes from amplifying effect of <span>\(Pr^{Beta}(r_i, r_j)\)</span>.</p>
<p>The idea of computation match result based on probability of winning one frame is not very novel. On <a rel="nofollow" href="https://wismuth.com/elo/calculator.html" target="_blank">this site</a>, authored by François Labelle, you can find online computation of probability of winning a “best of <span>\(N\)</span>” match (with some other functionality). I was very glad to notice that results of our computations match. However, I didn’t find any sources of incorporating this into Elo updating procedure. I will highly appreciate any information on this topic.</p>
<p>I’ve only managed to found <a rel="nofollow" href="http://freerangestats.info/blog/2015/08/07/fibs-elo-ratings-basics" target="_blank">this post</a> and <a rel="nofollow" href="http://www.bkgm.com/articles/McCool/ratings.html" target="_blank">this description</a> of Elo system on backgammon internet server (FIBS). Here different matches are handled by multiplying rating difference by square root of match length. However, there seems to be no strong theoretical reason for this.</p>
</div>

<div>
<h1>Experiment</h1>
<p>The experiment has several goals. Based on snooker data:</p>
<ul>
<li>Derive best value of “K factor” for both models.</li>
<li>Study stability of models in terms of prediction probability accuracy.</li>
<li>Study the effect of using “invitational” events on ratings.</li>
<li>Produce “best” rating history of all professional players from 2017/18 season.</li>
</ul>
<div>
<h2>Data</h2>
<details>

<div>
<pre><code># Function to split cases between &quot;train&quot;, &quot;validation&quot;, and &quot;test&quot; types
split_cases &lt;- function(n, props = c(0.5, 0.25, 0.25)) {
  breaks &lt;- n * cumsum(head(props, -1)) / sum(props)
  id_vec &lt;- findInterval(seq_len(n), breaks, left.open = TRUE) + 1
  
  c(&quot;train&quot;, &quot;validation&quot;, &quot;test&quot;)[id_vec]
}

pro_players &lt;- snooker_players %>% filter(status == &quot;pro&quot;)

# Matches only between pro players.
pro_matches_all &lt;- snooker_matches %>%
  # Using actually happened matches
  filter(!walkover1, !walkover2) %>%
  # Filter matches only between pro players
  semi_join(y = pro_players, by = c(player1Id = &quot;id&quot;)) %>%
  semi_join(y = pro_players, by = c(player2Id = &quot;id&quot;)) %>%
  # Add 'season' column
  left_join(
    y = snooker_events %>% select(id, season), by = c(eventId = &quot;id&quot;)
  ) %>%
  # Ensure arranging by end date
  arrange(endDate) %>%
  # Prepare for widecr format
  transmute(
    game = seq_len(n()),
    player1 = player1Id, score1, player2 = player2Id, score2,
    matchId = id, endDate, eventId, season,
    # Compute match type (&quot;train&quot;, &quot;validation&quot;, or &quot;test&quot;) with 50/25/25 split
    matchType = split_cases(n())
  ) %>%
  # Convert to widecr format
  as_widecr()

# Matches only between pro players in not invitational events (which by
  # quantity is dominated by Championship League).
pro_matches_off &lt;- pro_matches_all %>%
  anti_join(
    y = snooker_events %>% filter(type == &quot;Invitational&quot;),
    by = c(eventId = &quot;id&quot;)
  )

# Split confirmation
get_split &lt;- . %>% count(matchType) %>% mutate(share = n / sum(n))

# This should give 50/25/25 split (train/validation/test).
pro_matches_all %>% get_split()
## # A tibble: 3 x 3
##   matchType      n share
##   </code></pre></div></details></div></div></ul></div></div>