We present a theoretical analysis of the Elo rating system, a popular method for ranking skills of players in an online setting.
In particular, we study Elo under the Bradley–Terry–Luce model and, using techniques from Markov chain theory, show that Elo learns the model parameters at a rate competitive with the state of the art.
We apply our results to the problem of efficient tournament design and discuss a connection with the fastest-mixing Markov chain problem.

§ INTRODUCTION

The Elo rating system is a popular method for calculating the relative skills of players (or teams) in sports analytics and particularly chess <cit.>.
It is based on a simple zero-sum update rule: if player $i$ beats player $j$, then the rating of player $i$ increases proportionally to the model probability that $i$ would lose to $j$, while the rating of $j$ decreases by the same amount. This amount depends on the previously estimated difference in skills between $i$ and $j$.

Despite their widespread popularity, Elo rating systems still lack a rigorous theoretical understanding <cit.>. Here, we take a probabilistic approach and study Elo under the well-known Bradley–Terry–Luce model (BTL) <cit.>. In this model, the probability $p_{i,j}$ that $i$ wins against $j$

$w_i/(w_i + w_j)$,
where $w_k$ is the strength of player $k$.
In Elo, this is usually reparametrised via $w_k = \mathrm{e}^{\rho_k}$:

\[
\mathrm{e}^{\rho_i} / \rbr{\mathrm{e}^{\rho_i} + \mathrm{e}^{\rho_j}}
1/\rbr{1 + \mathrm{e}^{\rho_j - \rho_i}}
\sigma(\rho_i - \rho_j),
\]

where $\rho_k$ is the true rating of $k$
\(
\sigma(z) \cq 1/(1 + \exp{-z})
\)
for $z \in \mbr$
is the sigmoid function.
In this setting,
after observing $i$ beat $j$, the corresponding  ratings $\xx_i$ and $\xx_j$ are updated as

\[
\xx_i \leftarrow \xx_i + \eta \sigma(\xx_j - \xx_i)
\Qand
\xx_j \leftarrow \xx_j - \eta \sigma(\xx_j - \xx_i),
\]

where the step-size $\eta > 0$ is chosen by the modeller.
The size of the update depends exponentially on the difference in ratings:
beating a much lower rated opponent does not change the ratings much.

The goal of the Elo rating system is to estimate the true ratings of $n$ players by observing results of matches between pairs of players. It is, therefore, aiming to solve the problem of ranking from pairwise comparisons. Compared with most algorithms in the area <cit.>, however, Elo  benefits from three qualities that help explain its popularity in real-world applications: simplicity, interpretability and the ability to update a ranking of the players in an online fashion.

The knowledgeable reader might have noticed that Elo's update rule is actually based on the gradient of the BTL log-likelihood: Elo can be interpreted simply as stochastic gradient descent with fixed step-size. Rather than studying Elo from a convex optimisation angle, however, we take a more probabilistic point of view: we assume at each time $t$, players $i$ and $j$ are selected to play against one another with probability $q_{i,j}$. This allows us to interpret Elo as a Markov chain over $\mbr^n$ and deploy powerful tools to study its behaviour. On the other hand, this approach also presents us with some challenges: Elo is not a reversible Markov chain and, while it has a unique stationary distribution,
assuming a minor and natural condition on

it does not converge to it in total variation <cit.>.

§.§ Our Results

Our main contribution (Theorem <ref>) shows that Elo ratings, averaged over time, well-approximate the true ratings of the players with high probability.

In order to avoid the potential of unbounded ratings, which are unrealistic in practice, we consider a variant of Elo in which ratings are capped, whilst maintaining their zero-sum property.

We obtain rates of convergence with respect to the number of observed matches that are competitive against the state of the art in the BTL literature.
In contrast to most other algorithms for the BTL model studied in the literature, Elo learns the parameters of a BTL model in an online fashion.
These rates of convergence is remarkable since Elo was originally conceived as a simple ranking system for chess players.
Furthermore, our approach is very robust, and also applies to a parallel set-up in which multiple games are played concurrently.

We also discuss the problem of tournament design: we assume we are given a tournament (comparison) graph $G = ([n],E)$, and we would like to choose the match-up probabilities $(q_{i,j})_{\{i,j\} \in E}$ so that Elo's convergence rate is maximised. We highlight a connection between this problem and that of finding the fastest mixing Markov chain on $G$ <cit.>,
where the corresponding Markov chain on the graph is either in discrete or continuous time depending on whether we optimise number of games or parallel rounds, respectively.

As far as we know,

this connection has not been made formally before

in the sequential set-up, and both the analysis and optimisation of the parallel set-up are new.
In <ref>,
we also provide experimental results that showcase the usefulness of our strategy.

§.§ Related Work

Despite Elo ratings being the standard ranking system in many sports analytics communities <cit.>, there is a scarcity of work  analysing Elo from a probabilistic perspective. In particular, we are aware only of <cit.>, and the unpublished notes <cit.>
by the same author,
in which Elo is studied as a Markov process and convergence in distribution to a unique stationary distribution is proved.

If we consider Elo as a technique to estimate the parameters of the BTL model, there is a wealth of recent literature on the topic by the machine learning community <cit.>. In contrast to our setting, however, previous work typically considers an offline scenario, in which  the ratings of the players are computed after all the scheduled matches have taken place.
By standard concentration inequalities, this allows one to obtain a very good approximation of the probability a player wins against their neighbours in the comparison graph. The goal is then to deduce a global ranking of the players from such (very good) local information. In our setting, Elo ratings are dynamically updated before a good local approximation is achieved, which makes the analysis more challenging and requires the use of powerful, but delicate, concentration inequalities for Markov chains.

A comparison between the rate of convergence for Elo ratings obtained in our work and the rate of convergence of other algorithms for the BTL model is discussed in <ref>.

Finally, we mention that the connection between Markov chains and stochastic gradient descent with fixed step size has been studied before, e.g., in <cit.>.

§ CONVERGENCE RATES OF ELO RATINGS

In this section, we state our main result on the convergence rate of the time-averaged Elo ratings, discuss related work and outline the most important parts of the proof.

Throughout the paper,
“$f \lesssim g$” means “$f = O(g)$”,
“$f \ll g$” means “$f = o(g)$”
“whp” means “with probability $1 - O(1/n)$”;
also, later, $\EX{\pi}{\cdot}$ indicates that $X^0 \sim \pi$.

§.§ Our Results

We start with the explicit definition of our Markov chain.

Let $M \in \mbr$ and $n \ge 2$.
Let $\rho \in [-M, M]^n$ with $\sum_k \rho_k = 0$.

Let $q$ be a distribution on unordered pairs in $[n]$.
Let $\eta \in (0, \tfrac14)$.

A step of $\Elo_M(q, \rho; \eta)$ proceeds as follows.

* Suppose that the current vector of ratings is $x \in \mbr^n$.

* Choose unordered pair $\SET{I,J}$ to play according to $q$:

\[
\PR{}{ \SET{I,J} = \SET{i,j} }
\Quad{for all}
i,j \in [n].
\]

* Suppose that Player $I$ beats $J$.

Update ratings $x_I$ and $x_J$:

\[
\leftarrow
x_I + \eta \sigma(x_J - x_I);
\quad
\leftarrow
x_J - \eta \sigma(x_J - x_I).
\]

* Orthogonally project the full vector of ratings to $[-M, M]^n \cap \SET{ x' \in \mbr^n \mid \sum_k x'_k = 0 }$.

Let $X^t_k$ denote the rating of Player $k$ at time $t$, and $\pi$ the equilibrium distribution on $\mbr^n$.
Denote by

\[
\cq
\tfrac1t
%	\!
\sumt{s=T}[T+t-1]
%	\!
\Quad{for}
k \in [n]
\Quad{and}
t, T > 0
\]

the time-averaged ratings.
We typically start from $X^0 = (0,...,0)$.
The (deterministic) time $T$ is a burn-in phase, which allows the Elo ratings to get `near' the true skills, after which we start averaging.

The projection step ensures the Elo ratings do not become too large. It is required for our analysis, but not usually implemented in practice.
Indeed, our experiments, discussed in <ref>, suggest that, as long as the step-size $\eta$ is small enough, Elo ratings remain bounded.
Algorithms which estimate BTL parameters typically require some sort of projection <cit.> or regularisation <cit.>. A simple and efficient algorithm to realise the orthogonal projection is presented in the Appendix.

We show, for a suitable choice of parameters $t$, $T$ and $\eta$, that the time-averages are concentrated around the true ratings.
In other words, we can use Elo to obtain an MCMC estimate of the true ratings.

Elo ratings in equilibrium are, in general, a biased estimator of the true ratings:
i.e., $\EX{\pi}{X^0_k} \ne \rho_k$.

Hence, there will be both a bias and an error term in our MCMC-type estimate of $\norm{ A^{t,T} - \rho }_2$.

The MCMC convergence rate depends on a spectral gap $\lambda_q$.
This parameter quantifies how fast local  information about the relative strengths of two players is propagated to the rest of the ratings. Similar parameters appear in most of the related literature, e.g., <cit.>.

Let $q$ be a distribution on unordered pairs in $[n]$.
$q_{i,j} \cq q_{\set{i,j}}$,
$d_{i,j} \cq \one{i = j} \sum_k q_{i,k}$
$\Delta_{i,j} \cq d_{i,j} - q_{i,j}$
for $i,j \in [n]$.
Let $\lambda_q$ denote the spectral gap
—i.e., the second smallest eigenvalue—
of the Laplacian $\Delta$.
Always, $\Delta \mathbf 1 = \mathbf 0$ and $\Delta \succeq 0$.

Equivalently, $\lambda_q$ is the spectral gap of the continuous-time Markov chain on $[n]$ with
rates $q_{i,j} = q_{\set{i,j}}$ for $i,j \in [n]$.
Note the scaling: $\sum_{i,j} q_{i,j} = 2$.
So, the typical time until the continuous-time chain jumps is order $n$, not order $1$.
This implies $\lambda_q \le 4/n$ always.

We assume $\lambda_q > 0$. This holds unless there exists a non-empty subset $S \subsetneq [n]$ with $\sum_{i \in S, j \not \in S} q_{i,j} = 0$—i.e., players in $S$ never play those in $S^c$.
This makes estimation of the ratings impossible.

Our main result measures the disparity between the time-averaged Elo ratings and the true ratings.

Let $X \sim \Elo_M(q, \rho; \eta)$, as in <ref>.
Let $C_1, C_2 < \infty$.
Then, there exists a constant $C_0$, depending only on $(C_1, C_2)$, such that if

\[
\min\SET{\lambda_q, \eta, 1/t} \ge n^{-C_1}
\Quad{and}
\min\SET{t, T} \ge C_0 \tmix,
\Qwhere
\tmix
\cq
e^{2M} \eta^{-1} \lambda_q^{-1} \log n,
\]

\[
\PR[\bigg]{}{
\tfrac1n \norm{ A^{t,T} - \rho }_2^2
\le
\frac{C_0 e^{4M}}{\lambda_q n}
\biggl( \eta + \frac{(\log n)^2}{\lambda_q t} \biggr)
\ge
1 - n^{-C_2}.
\]

In particular, if $\eta \asymp (\log n)^2 / (\lambda_q t)$, then

\[
\tfrac1n \norm{ A^{t,T} - \rho }_2^2
\lesssim
\frac{e^{4M} (\log n)^2}{\lambda_q n}
\frac1{\lambda_q t}
\quad
\whp
\Quad{as}
n \to \infty.
\]

Ideally, $\lambda_q n = \tilde \Omega(1)$;
e.g., if $q$ is uniform over the edges of an expander graph, then $\lambda_q n \asymp 1$.
In this case, we view $(\log n)^2 / (\lambda_q n) = \tilde O(1)$

as the error's pre-factor and $1/(\lambda_q t)$ as the (squared) convergence rate.

Also, we can then choose $\eta = \tilde \Omega(1)$ and obtain non-trivial convergence results.

This choice of $\eta$ is comparable to that used in practice.

on average, only $\tilde O(1)$ games per player
are need to be observed to guarantee good approximation of the ratings.

The $\eta$ term arises from the average bias $\tfrac1n \norm{ \EX{X\sim\pi}{X} - \rho }_2^2$, which is non-zero in general.
An estimate on the average variance $\tfrac1n \sum_k \VAR{X\sim\pi}{X_k}$ is necessary to obtain the MCMC convergence rate.

If $\pi$ is the equilibrium distribution of $\Elo_M(q, \rho; \eta)$,

\[
\tfrac1n
\norm{ \EX{X\sim\pi}{X} - \rho }_2^2
\le
4 e^{4M} \eta / (\lambda_q n)
\Qand
\tfrac1n
\sumt{k} \VAR{X\sim\pi}{X_k}
3 e^{2M} \eta / (\lambda_q n).
\]

Despite Elo ratings being biased, the estimated probability a player wins their next match is not:

\[
\sumt{j \in [n]}
q_{i,j} \EX{X\sim\pi}{ \sigma(X_i - X_j)}
\sumt{j \in [n]}
q_{i,j} \sigma(\rho_i - \rho_j)
\Qforall
i \in [n].
\]

This holds when no projection step is performed as part of the Elo update (equivalently, $M \cq +\infty$).

The mixing time is more delicate, as the chain makes deterministic-size discrete jumps, but its equilibrium distribution is continuous and, therefore, the chain does not converge in total variation.
Instead, we measure convergence in the Wasserstein, also known as transportation, distance.

If $X, Y \sim \Elo_M(q, \rho; \eta)$, then there exists a step-by-step coupling with

\[
\EX{(x,y)}{ \norm{X^t - Y^t}_2^2 }
\le
(1 - \kappa)^t
%	\rbb{ 1 - \tfrac1{16} e^{-2M} \eta \lambda_q }{}^t
\norm{x - y}_2^2
\Qwhere
\kappa \cq \tfrac18 e^{-2M} \eta \lambda.
\]

Markov chains satisfying the contraction property (i.e., positively curved) satisfy powerful concentration inequalities,
developed particularly by Joulin and Ollivier <cit.>.

The idea is that if $\norm{X^0 - Y^0}_2 = D^0$, then $D^t \cq \norm{X^t - Y^t}_2$ is small when $t \asymp \kappa^{-1} \log D^0$.

For us, $D^0 \le 2 M n$, leading to $\tmix \asymp \kappa^{-1} \log n$.
In the reversible case, this can provide bounds, e.g., on the spectral gap.

The approach is particularly applicable to Markov chains on finite metric spaces, such as graphs.
The set-up of finite graphs was analysed by Bubley and Dyer <cit.> under the name path coupling.
In this case, $\ex{D^t} \le \tfrac12$ implies $\pr{X^t = Y^t} \ge \tfrac12$, since the minimum graph distance between $x \ne y$ is $1$.

Our underlying metric space is $(\mbr^n, \norm\cdot_2)$, however.
Inequalities for this more general set-up have been developed, but they often give weaker bounds.
In our case, they lose a factor $1/(\eta \sg) \gtrsim n$ in the convergence rate.
Morally, though, the exponential convergence at rate $\kappa$ still implies that $1/\kappa$ is the correct timescale for mixing and concentration,
perhaps up to some logarithmic factors.
Establishing this rigorously is the most challenging part of our work from a technical point of view.

§.§ Comparison with Related Work

In this section, we compare the convergence rate of Elo given by <ref> against the state of the art for BTL estimation.
We highlight, once again, that previous work focusses on the problem of offline estimation, whilst Elo works in the more challenging online setting.
Nevertheless, Elo is able to match the state-of-the-art algorithms in the offline setting for a wide range of parameters.

Our results imply that Elo provides an estimator $\hat{\rho}$

of $\rho$
$\|\hat{\rho} - \rho\|_2^2 \lesssim \mathrm{e}^{4M} (\log n)^2 / (\lambda_q^2  t)$ whp.
This matches, up to a log factor, the results by
Hajek, Oh and Xu
<cit.>, who prove that the MLE constrained on
$[-M, M]^n \cap \{x \in \mbr^n \mid \sum_k x_k = 0\}$,

$\rho_{\text{MLC}}$, satisfies
$\|{\rho}_{\text{MLC}} - \rho\|_2^2 \lesssim \mathrm{e}^{8M} \log n / (\lambda_q^2 t)$.
The dependency on $\lambda_q$ has been improved by
Shah, Balakrishnan and Bradley
<cit.>, who show $\|{\rho}_{\text{MLC}} - \rho\|_2^2 \lesssim n \mathrm{e}^{8M} \log n / (\lambda_q t)$. Since $n/\lambda_q \ge 1$, this is at least as good as our result and potentially better for certain choices of $q$. However, up to log factors, our result matches theirs when $q$ corresponds to sampling edges of an expander graph—e.g., a complete graph, or an Erdős–Rényi graph with parameter $p \gg (\log n)/n$.
Our result improves

the constant in front of $M$.

More recently,
Li, Shrotriya and Rinaldo
<cit.> proved a regularised version of the MLE, $\rho_{\text{MLR}}$, achieves $\|{\rho}_{\text{MLR}} - \rho\|_2^2 \lesssim \mathrm{e}^{2M_E} \delta / (\lambda_q^2 t)$, where $\delta$ is the ratio between the maximum and average degree of the comparison graph $G=([n],E)$ with $E = \{\{i,j\} \mid q_{i,j} > 0\}$ and $M_E = \max_{i,j : q_{i,j} > 0} \abs{\rho_i - \rho_j} \le 2 \max_k \abs{\rho_k}$.
A version of our analysis also applies with $M_E$, instead of $M$, but we are not able to prove $\max_{i,j : q_{i,j} > 0} \abs{X^t_i - X^t_j} < M_E + 1$ holds for polynomially long.

Notice that their result is weaker when the comparison graph is relatively sparse but has a few high degree nodes.
Moreover, they require each player to play the same number of matches;
i.e., $q_{i,j} = 1/\abs E$ for $\set{i,j} \in E$.
Both limitations are particularly problematic in the context of tournament design discussed in <ref>.
On the other hand, they obtain $\ell_\infty$ error bounds too.
We refer to <cit.> for
further discussion of previous work.

§.§ Outline of Proof

We now highlight the key steps in estimating the MCMC-type error $\tfrac1n \norm{ A^{t,T} - \rho }_2^2$ in <ref>, where $A^{t,T}$ is the $t$-step time average of the ratings $X$ after a burn-in phase of length $T$.
We use

\[
\norm{ A^{t,T} - \rho }_2
\le
\norm{ A^{t,T} - \ex[\pi]{X^0} }_2
+	\norm{ \ex[\pi]{X^0} - \rho }_2,
\]

bounding these error and bias terms separately,
via Theorems <ref> and <ref>, respectively.

Their proofs rely on estimating the change in $\ell_2$ norm a single step, and using the Lipschitz property of the sigmoid function $\sigma$, along with some other careful manipulations.
The full proofs are given in the appendix.

The primary challenge is to leverage the positive curvature result of <ref> to deduce the required concentration inequality in <ref>.
As noted at the end of <ref>, positive curvature is well-suited to this goal, but the specifics of our set-up cause significant difficulties.
Particularly, the aforementioned results reliant only on curvature are not sufficiently strong for our purposes.

A further complication is that Elo is a non-reversible Markov chain.

There is a general understanding, or perhaps belief, in the MCMC community that non-reversibility can improve concentration results.
This has lead to many strategies being proposed,
such as
see, e.g., <cit.>.
Even an alternative definition of the spectral gap has been proposed <cit.>.

Unfortunately, the majority of general concentration results based on spectral properties <cit.> are actually worse in the non-reversible set-up:
they bound the convergence rates by those of the multiplicative or additive reversibilisation.

These can be hard to estimate, needing detailed knowledge of the equilibrium distribution, and the resulting bounds are often crude.

Finally, we mention that bounds on the total-variation mixing time can be leveraged to provide concentration bounds; see, particularly, <cit.>.

However, Elo does not converge in total variation!
Indeed, $X^t$ is finitely supported, on at most $(2n)^t$ states, given $X^0$, but $\pi$ is continuous.
These $(2n)^t$ states can be used to distinguish $X^t$ and $\pi$.
Nevertheless, it is this mixing-time approach that we adapt.

The contraction property implies that two realisations $X$ and $Y$ can be coupled so that their expected relative distance decreases exponentially, with rate $\kappa$.
This is ideal for using the coupling approach to mixing.
However, it is not possible to guarantee that $X^t = Y^t$ at some point, due to the chain's finite support.
Instead, we analyse a noisy version:
after each update, some (continuous) noise is added to the ratings.
This must be added carefully to preserve the contraction in <ref>.

Denote the noisy versions $U$ and $V$.
Once $U$ and $V$ are sufficiently close, the additive noise can be used to couple them exactly, to get $U^t = V^t$.
The size of the noise must be balanced:
too small and this final coupling step is too difficult;
too big and $U$ can no longer be compared with $X$.

We delve deeper into this noisy approximation, explaining what noise to add, how it accumulates and what the resulting mixing time is.
The interaction of the noise with the projection step is delicate and technical; we omit it from the description here, but it should be kept in the back of the mind.

We consider the following noisy version $U$ of the Elo process $X$.
Let $\delta > 0$.

Draw $U^{t+1/2}$ according to an uncapped Elo step started from $U^t$
—i.e., as if $M = \infty$.
Suppose that i j were chosen
—i.e., $\set{ k \mid U^{t+1/2}_k \ne U^t_k } = \set{i,j}$.

* Draw $\tilde U_i, \tilde U_j \sim^\iid \Unif([- \sqrt \delta, + \sqrt \delta])$
independently and set $\tilde U_k \cq 0$ for $k \notin \set{i,j}$.

\[
\cq
U^{t+1/2}_k + \tilde U_k
\Qfor
k \in [n].
\]

This does not preserve the zero-sum property, as the additive noise is independent.
The independence is crucial later to allow two version to coalesce.
We need to control the cumulation until that point.

We control the difference between $X$ and $U$ via the natural coupling:
pick the same pair of players to play, and observe the same result;
sample the noise independently.

Careful computation gives

\[
\norm{ X^{t+1} - U^{t+1} }_1
\le
\norm{ X^{t+1/2} - U^{t+1/2} }_1
+	2 \sqrt \delta
\le
\norm{ X^t - U^t }_1
+	2 \sqrt \delta
\le
\le
2 t \sqrt \delta.
\]

The second inequality actually says that Elo is non-negatively curved in $\ell_1$.
Iterating this,

\[
\normb{
\tfrac1t \sumt{s=0}[t-1] X^s
-	\tfrac1t \sumt{s=0}[t-1] U^s
\le
\tfrac1t
\sumt{s=0}[t-1]
\norm{ X^s - U^s }_1
\le
t \sqrt \delta
\quad
\text{deterministically}.
\qedhere
\]

We use the natural coupling between two noisy versions $U$ and $V$:
choose the same players,
observe the same result
add the same noise.
Then, there is still rate-$\kappa$ contraction:

\[
\ex[u^0,v^0]{ \norm{ U^1 - V^1 }_2 }
\le
\ex[u^0,v^0]{ \norm{ U^{1/2} - V^{1/2} }_2 }
\le
(1 - \kappa)
\norm{ u^0 - v^0 }_2.
\qedhere
\]

We bound the total-variation mixing time of the noisy chain via the coupling method:

\[
\normb{ \pr[u^0]{U^t \in \cdot} - \pr[v^0]{V^t \in \cdot} }_\textup{TV}
\le
\pr[(u^0, v^0)]{ U^t \ne V^t }.
\]

First, we burn-in using the natural coupling,
then use a new coupling which exploits the noise.

We start with the natural coupling as given above.
Using the rate-$\kappa$ curvature,

\[
\ex{ \norm{ U^t - V^t }_2 }
\le
(1 - \kappa)^t
\norm{ U^0 - V^0 }_2
\le
2 M n e^{-\kappa t}
\le
\delta^2
\Qif
\ge
\cq
\kappa^{-1} \log(2Mn / \delta^2).
\]

\[
\pr{ \norm{ U^t - V^t }_\infty > \delta }
\le
\pr{ \norm{ U^t - V^t }_2 > \delta }
\le
\delta
\Qif
\ge
\]

Once the absolute difference $\abs{ U^{t_\delta}_k - V^{t_\delta}_k }$ are at most $\delta$ for all $k$,
the additive noise, which is order $\sqrt \delta$,
dominates the change in difference after a single Elo step, which is order $\delta$.
Thus, with complementary probability order $\delta / \sqrt \delta = \sqrt \delta$, we can couple the noise so that the rating of a player is the same in both $U$ and in $V$ after they play a game.
Moreover, such a successful step preserves the $\ell_\infty$ bound of $\delta$.

The $\ell_\infty$ bound implies that all players are chosen within $t_\delta$ steps with probability at least $1 - \delta$.
Hence, the probability of not successfully matching all ratings after $t_\delta$ steps is at most order $\delta + \sqrt \delta t_\delta$.

Combining these two bounds gives a bound of order $\delta + \sqrt \delta t_\delta$ on the total-variation distance at time $2 t_\delta$.
The polynomial-growth/decay assumptions in the theorem allow us to choose $\delta$ to be an appropriate inverse polynomial (in $n$) and obtain $\delta + \sqrt \delta t_\delta \ll 1$ and $t_\delta \asymp \kappa^{-1} \log n \asymp \tmix$.

The above mixing-time bound allows us to establish concentration of time-averages of the noisy Elo ratings, using results from <cit.>.
Particularly, under certain assumptions, if $Z$ is a Markov chain with equilibrium distribution $\pi$ and mixing time $\tmix$, then

\[
- \log \pr*[\pi]{ \absb{ \tfrac1t \sumt{s=0}[t-1] f(Z^s) - \pi_f } \ge \zeta }
\gtrsim
\sigma_f^{-2} \zeta^2 t / \tmix,
\]

where $\pi_f = \ex[\pi]{f}$ and $\sigma_f^2 \cq \var[\pi]{f}$.
Notice that this requires $Z^0 \sim \pi$, which we do not impose; this requirement can be circumnavigated by using a burn-in, again comparing with

a noisy version.

We want to take $f \cq f_k$ to be the projection onto the $k$-th coordinate, which is the $k$-th player's rating, for each $k$, then do a union bound over the $n$ players.
This motivates taking

\[
\zeta_k^2
\asymp
\sigma_{f_k}^2 \tmix
\frac{\log n}{t},
\Quad{which satisfies}
\tfrac1n
\sumt{k}
\zeta_k^2
\asymp
\frac{\eta}{\sg n}
\frac{\log n}{\sg \eta}
\frac{\log n}{t}
\frac1{\sg n}
\frac{(\log n)^2}{\sg t},
\]

using <ref>.
This is the decay rate required in <ref>,
but for the noisy version.

We need to compare the noisy version with the original.
We do this via the estimate established earlier:

\[
\normb{
\tfrac1t \sumt{s=0}[t-1] X^s
-	\tfrac1t \sumt{s=0}[t-1] U^s
\le
\tfrac1t
\sumt{s=0}[t-1]
\norm{ X^s - U^s }_1
\le
t \sqrt \delta
\quad
\text{deterministically}.
\]

We are taking $t$ to be polynomial in $n$, and can choose $\delta^{-1}$ to be a sufficiently large polynomial so that this accumulated error is small.
Care must be taken to to make all this rigorous, but once it is done, we are able to deduce the concentration result for the original Elo process.

§ TOURNAMENT DESIGN

§.§ Our Results

In this section, we assume we are given a comparison graph $G = ([n], E)$, where edge $\SET{i,j} \in E$  indicates that Players $i$ and $j$ are able to play against one another.
We want to construct a distribution $q$ over edges $E$ of $G$ so that

Elo can most efficiently approximate the true ratings of the players.

The decay rate in <ref> is governed by the spectral gap $\sg$.
So, we want to maximise $\sg$.

\[
\textstyle
\lambda^\star_\cts
%	\lambda^\star_\cts(G)
\cq
\sup\SET[\big]{ \lambda_q \given q \in [0,1]^E, \: \sum_{e \in E} q_e = 1 }.
\]

This equals the largest spectral gap achievable by a continuous-time Markov chain on $[n]$ with
transitions only across edges of $G$
average jump-rate $1/n$.
It is the fastest-mixing Markov chain problem, introduced by Sun et al. <cit.>, and can be formulated as a semidefinite program.
Olesker-Taylor and Zanetti <cit.> recently proved that

\(
%	1/(\lambda_q n)
\lambda^\star_\cts n
\gtrsim
1 / (\operatorname{diam} G)^2.
\)
This implies the following.

Suppose that the comparison graph $G = ([n], E)$ is given.
In the set-up of <ref>,
there exists a distribution $q$ on the edges $E$ such that

\[
\tfrac1n
\norm{ A^{t,T} - \rho }_2^2
\lesssim
n (\operatorname{diam} G)^2 (\log n)^2 / t
\quad
\whp
\Qif
\eta
\lesssim
n (\operatorname{diam} G)^2 (\log n)^2 / t.
\]

This choice of probabilities $q$ can improve drastically over uniform weights—as used by <cit.>.
E.g., if $G$ consists of two cliques connected by $k$ edges, then $\lambda_q n \lesssim k/n^2$ when $q_e \cq 1/\abs E$ is uniform over $E$, whilst the optimal $\lambda^\star_\cts n \asymp 1$ <cit.>.

As a consequence, with the optimal choice of $q$, only $\tilde{O}(n)$ total matches need to be played to obtain a good approximation of the true ratings, compared with $\tilde{O}(n^3/k)$ for a uniform $q$.

This demonstrates the power of being able to choose $q$, given $G$.

Elo naturally parallelises:
if two games consist of disjoint pairs of players,

then the Elo update resulting from one is independent of the result of the other.
Hence, if we wish to minimise the number of rounds (i.e., sets of games that can be played in parallel), rather than games,
we should consider a distribution $\tilde q$ on the set $\mcm$ of matchings of the graph $G = ([n], E)$—i.e., collections of disjoint edges.

Let $\tilde q$ be a distribution on $\mcm$.
A single step of $\ParElo_M(\tilde q, \rho; \eta)$

first selects a matching $S \subseteq E$ according to $\tilde q$ and applies the Elo update (with scale $\eta$) to each pair in $S$.
Then, the resulting vector is orthogonally projected back to zero-sum vectors in $[-M,+M]^n$.

Our analysis is robust enough to handle the parallel case, obtaining convergence rate $1/(\lambda_q t)$, where now $q_e$ is the marginal probability that edge $e \in E$ appears in the matching:

\[
\textstyle
\cq
\sum_{S \in \mcm : e \in S}
\tilde q_S.
\]

Let $X \sim \ParElo_M(\tilde q, \rho; \eta)$.

Then, under the conditions of <ref>,

\[
\tfrac1n \norm{ A^{t,T} - \rho }_2^2
\lesssim
\frac{e^{4M} (\log n)^2}{\lambda_q n / N}
\frac1{\lambda_q t}
\quad
\whp
\Qif
\eta
\asymp
\frac{(\log n)^2}{\lambda_q t},
\]

$N \cq \sum_{e \in E} q_e$ is the mean size of the matching.
To emphasise, here, $t$ and $T$ count rounds.

It can be shown that this factor $N$ is needed in the pre-factor via a time-change analysis.

It is natural to optimise $\lambda_q$ over $q$ which can arise as the marginals of a distribution $\tilde q$ on $\mcm$.

\[
\cq
\sumt{e \in E : k \in e}
\le
\Quad{for all}
k \in [n]
\]

for any such $q$.
Thus, the matrix $Q = (q_{i,j})_{i,j \in V}$

is substochastic, so $\lambda_q$ corresponds to the spectral gap of the discrete-time Markov chain on $G = ([n], E)$ with weights $q_{i,j} = q_{\SET{i,j}}$.

\[
\textstyle
\lambda^\star_\disc
\cq
\sup\set{ \lambda_q \given q \in [0,1]^E, \: \max_{k \in [n]} q_k \le 1 }.
\]

This is the optimal spectral gap of a discrete-time Markov chain on $[n]$ with
transitions allowed only across edges of $G$
uniform stationary distribution.
Again, $\lambda^\star_\disc$ can be formulated as a semidefinite program <cit.>
and is related to the vertex conductance via a Cheeger-type inequality <cit.>.

We show that a distribution $\tilde q$ over matchings with $\lambda_q \ge \tfrac13 \lambda^\star_\disc$ can be found by decomposing the substochastic $Q$ into a convex combination of permutation matrices using the Birkhoff–von-Neumann theorem, followed by decomposing each permutation into disjoint cycles.

Suppose that the comparison graph $G = ([n], E)$ is given.
In the set-up of <ref>,
there exists a distribution $\tilde q$ on matchings $\mcm$ such that

\[
\tfrac1n \norm{ A^{t,T} - \rho }_2^2
\lesssim
\frac{e^{4M} (\log n)^2}{\lambda^\star_\disc n / N}
\frac1{\lambda^\star_\disc t}
\quad
\whp
\Qif
\eta
\asymp
\frac{(\log n)^2}{\lambda^\star_\disc t}.
\]

In many examples, such as if $G = ([n], E)$ is an expander, but also if $G$ is a cycle, then $\lambda^\star_\disc \asymp n \lambda^\star_\cts$.
In this case, parallel Elo really is as good, up to constants, as $n$ steps of the original (`series') Elo.
(Recall that $N = \sum_e q_e = 1$ is required for $\lambda^\star_\cts$, but that $N \asymp n$ is possible for $\lambda^\star_\disc$.)

Other times, there is already a continuous-time chain,
with average jump-rate $1$,
which has spectral gap order $1$;
e.g., two cliques connected by a single edge.
In this case, the optimal parallel version gives no real improvement,
even measured by rounds,

over the optimally weighted series version.

§.§ Comparison with Related Work

The problem of designing an efficient tournament graph is related to active ranking <cit.>. Active ranking, however, allows one to choose which matches to schedule next after observing the results of some matches.

In contrast,
we design a probability distribution over matches offline,
without the possibility of changing such distribution after observing some results.

This problem has been considered by
Li, Shrotriya and Rinaldo
They discuss a divide-and-conquer strategy that essentially requires oversampling edges across bottlenecks in the graph.
The drawback of this strategy is that it requires partitioning the graph into well-connected pieces, which is a non-trivial task itself.
Furthermore, <cit.> does not provide explicit bounds on the sample complexity that can be obtained in this way, besides discussing a few examples where the bottleneck is known.

Nonetheless, our approach shares some similarity with <cit.>: the fastest-mixing Markov chain implicitly up-weights edges across bottlenecks.
The main advantage, however, is that it provides the optimal spectral gap allowed by a graph topology, which is related to the diameter of the graph <cit.>.
Moreover, both $\lambda^\star_\cts$ and $\lambda^\star_\disc$ can be formulated as SDPs, for which fast (polynomial-time) solvers exist.
We remark, however, that this construction might require certain nodes to play an overwhelmingly large number of matches: this would be problematic for the results of <cit.>.
As far as we know, minimising the number of `parallel rounds' has not been considered before.

§ EXPERIMENTAL RESULTS

We close with discussion of some specific examples and experimental results.
Additional experiments are discussed in the Appendix.

We start by considering dumbbell tournament graphs consisting of two cliques of $n/2=20$ vertices connected by a matching of $k \in \{1,20\}$ edges.
For each graph, we perform Elo simulations where the match-ups between players are sampled according to the following probability distributions.

* The uniform distribution $q_{\mathcal{U}}$ over the edges of the graph.

* The optimal sequential $q_{\text{seq}}^\star$ derived from the fastest-mixing continuous-time Markov chain.

* The distribution over matchings $q_{\text{par}}^\star$ derived from the fastest mixing discrete-time Markov chain,
where multiple games are played in parallel in each round.

From <ref>,
we know that
$\lambda_{q_{\mathcal{U}}} \asymp k/n^3$, $\lambda_{q_{\text{seq}}^\star} \asymp 1/n$ and $\lambda_{q_{\text{par}}^\star} \asymp k/n$.

We sample the true ratings of the players according to independent Gaussians, with mean equal to 1 on one clique, mean 2 on the other, and standard deviation equal to $0.2$ in both. This difference in average ratings between cliques simulates a scenario where the two cliques correspond to two different leagues of slightly different strength on average.

We perform Elo simulations initialising the Elo ratings at zero and setting $\eta = 0.1$. For each graph, we repeat each experiment ten times, each time sampling new true ratings. Experimental results are displayed in Figures <ref> and  <ref>. Simulations are repeated ten times: lines correspond to the average $\ell_2^2$-distance between time-averaged Elo ratings and the true ratings, divided by the number of players ($n = 40$). Shaded regions corresponds to 25–75 percentile over these ten trials. Cyan and blue lines correspond to the same experiments where we have sampled multiple games in parallel as described in <ref>;
we display the decay of error wrt the total number of games played in cyan and wrt the number of parallel rounds in blue.

The blue line is solid for the first $2\cdot 10^4$ games (same number of games displayed in cyan).

We display the experimental results for $k=1$, i.e., two cliques of $20$ vertices connected by a single edge, in Figure <ref>. The  experiments align with the theoretical results of  <ref>: when $k=1$, the spectral gap $\lambda_{q_{\text{seq}}^\star}$ of the optimal sequential distribution is of the same order of the spectral gap for our nearly parallel construction $\lambda_{q_{\text{par}}^\star}$. Indeed, if convergence is measured wrt the number of rounds, the two corresponding errors seems to decay at a similar same rate, with the parallel version being slightly better, but requiring more than ten times the total number of games. As predicted by the theory, both distributions result in much faster convergence than the uniform distribution. If we measure the total number of games rather than rounds, the parallel version still results in a faster convergence than the uniform distribution, but not overwhelmingly so.

[width=,,trim=.02 0 .09 .05,clip]figs/dumbbell_n40_b=1.png
Elo simulation results for a dumbbell graph with one edge between two cliques of $20$ vertices. Match-ups are sampled from three different probability distributions.

[width=,,trim=.02 0 .09 .05,clip]figs/dumbbell_n40_b=20.png
Elo simulation results for a dumbbell graph with a perfect matching of $20$ edges between two cliques of $20$ vertices.

We now consider a dumbbell graph with $k=n/2$, i.e., two cliques connected by a perfect matching. In this case, the fastest discrete- and continuous-time Markov chains have the same order-$1$ spectral gap. However, since the spectral gap for sequential Elo is then rescaled by a factor of $1/n$, to achieve convergence, we expect the number of rounds for the parallel version to be much smaller than the number of games required by the sequential one; convergence should instead happen at the same rate when measured in the total number of games. This is clearly shown in Figure <ref>. The uniform distribution performs much worse, which we would expect from its order-$1/n^2$ spectral gap.
In the optimal parallel and sequential distributions, the probability to sample an edge from the bottleneck or from inside the cliques is balanced, while the uniform distribution oversamples edges inside the cliques.
Experiments for the intermediate case of $k \in \{5,10\}$ are discussed in the Appendix.

[width=,,trim=0 0 0 0,clip]figs/max_complete.png
Largest Elo rating in absolute value for a complete graph of varying size. True ratings are uniformly distributed in $[-1,1]$.

We end this section by discussing experimental results concerning the maximum rating reached by Elo. Figure <ref> displays the behaviour of the largest Elo rating in absolute value, where pairs of players are selected uniformly at random (i.e., the underlying graph is a complete graph). The true ratings are sampled uniformly at random in $[-1,1]$.

The initial Elo ratings are set equal to zero. We simulate up to 50000 matches for a number of players that goes from 100 to 1000.
We observe that the maximum rating is always below 1.75, corroborating our belief that, in many scenarios,
the maximum Elo rating
does not diverge for a long time.
Further experiments and discussions are given in the Appendix.

§ CONCLUSION AND OPEN PROBLEMS

Our work is a first step towards establishing
the theoretical foundations of
the Elo rating system, a popular ranking method in sports analytics.
In particular, our main contribution is an analysis of Elo under the BTL model, establishing convergence results competitive with the state of the art.

There are several questions prompted by our work.
from a technical point of view,
we would like to control the maximal Elo rating.

This is necessary to understand when we can remove the projection step in our definition of the Elo Markov chain, and make our theoretical results more aligned with practice.

Moreover, we would like to better understand the shape of the stationary distribution
of Elo.
This could help us, for example,
obtain bounds on the rate of convergence in $\ell_\infty$.

Finally, a touted strength of the Elo rating system in practical applications is its ability to dynamically update the ratings in response to changes in players' skills.
Can we model these changes in a way that allows us to prove Elo can keep track of them, and bound the corresponding mean squared error?
Such questions are often studied in the statistics literature;
see, e.g., the recent paper <cit.> for details.

§ EXPERIMENTS: FURTHER DISCUSSION AND FIGURES

In this appendix we discuss additional experimental results. Code for our experiments is included in the supplementary material.

In <ref> we display further experiments related to dumbbell graphs.
In particular, we consider graphs consisting of two cliques of $20$ vertices connected by $5$ (left)
and $10$ (right) edges arranged in a matching.
The experimental setup is the same as the one discussed in <ref>. Notice how, the more edges are added to the bottleneck,
the better the distribution optimised for parallel rounds performs. The optimal sequential distribution, instead, performs roughly the same:
adding edges in the bottleneck doesn't really improve its spectral gap, which depends mainly on the (unchanged) diameter.

[width=0.4,,trim=.02 0 .09 .05,clip]figs/dumbbell_n40_b=5.png
[width=0.4,,trim=.02 0 .09 .05,clip]figs/dumbbell_n40_b=10.png

simulation results for dumbbell graphs with $k=5$ (left) and $k=10$ (right) edges between two cliques of $20$ vertices.

We also consider a pyramidal graph, which is constructed as follows. We first sample three Erdős–Rényi random graphs with size, resp., $n_1 = 64$, $n_2=32$ and $n_3 = 16$, and density $p = 1/2$. We then connect the first graph to the second and the second to the third with two sparse cuts (see Figure <ref>). This graph is constructed to loosely resemble the pyramidal structure of, e.g., national sport leagues.

Schematic representation of the pyramidal graph.

[width=,,trim=.02 0 .09 .05,clip]figs/pyramid.png
Elo simulation results for the pyramidal graph.

We conduct Elo simulations with the same set-up as for dumbbell graphs discussed earlier. In particular, we repeat the simulations ten times, resampling each time the players' true ratings. The true ratings are sampled as follows: independent normal distributions of standard deviation 0.2 and mean $0$ for the Erdős–Rényi at the bottom of the pyramid, mean $1$ for the Erdős–Rényi in the middle, and mean $2$ for the one at the top. This is, again, to loosely simulate the fact that sport leagues are characterised by stronger players/teams towards the top of the pyramid. Experimental results are shown in Figure <ref>. In particular, we observe that the rate of converge for the uniform distribution is much slower than for the optimal sequential one. This is, again, predicted by the results of <ref>: the two bottlenecks slow down the convergence in the uniform case; while the small diameter assures faster convergence for the optimal distribution by Corollary <ref>. The distribution optimised for parallel Elo, instead, guarantees very fast convergence when measured according to the number of rounds.

In <ref> we display experimental results for a topology corresponding to the giant component of an Erdős–Rényi random graph of $n=100$ vertices and density $p=0.02$. True ratings are distributed as independent standard Gaussians. This graph has reasonably good connectivity, so we expect the number of games required to converge to be roughly equal for all the sampling distributions considered, with a good scope for parallelisation. This is confirmed by the error plot.
What is perhaps surprising is that the optimal sequential distribution actually performs slightly worse than the uniform one. This can be explained by the fact that the spectral gap measures the rate of convergence for the worst possible vector of ratings: if the ratings do not depend on the topology of the graph, like in this example, the inverse of the spectral gap is an overtly pessimistic upper bound. Indeed, the optimal sequential distribution tends to oversample nodes that are in a central position in the graph and undersample nodes at the periphery. This should make convergence faster because it allows faster movement of information between distant nodes, but in this specific example is not helpful: since ratings are sampled in a iid fashion, the
faster movement of information globally is not very helpful and is upset by slower convergence for undersampled nodes.

[width=0.45,,trim=0 0 0 0,clip]figs/ER_n=100_p=0.02_topo.png
[width=0.4,,trim=.02 0 .09 .05,clip]figs/ER_n=100_p=0.02.png

Left: topology of the giant component of an Erdős–Rényi random graph of $n=100$ vertices and density $p=0.02$; edges are reweighed according to the distribution corresponding to the fastest mixing continuous-time Markov chain.
Right:  simulation results for the same graph.

We now present further experiments about the maximum Elo rating (in absolute value). Again, we sample true ratings independently and uniformly in $[-1,1]$. In <ref>, we present experiments for a path (left) and star (right) topology of varying size. In both cases, the largest Elo rating in absolute value remains smaller than twice the largest true rating. Notice that in the star graph the central node plays all the matches; this means the same player plays in total fifty-thousands matches without its value becoming particularly large.

[width=0.4,,trim=0 0 0 0,clip]figs/max_path.png
[width=0.4,,trim=0 0 0 0,clip]figs/max_star.png

Behaviour of the largest Elo rating in absolute value for path (left) and star (right) graphs of varying size.

Of course, these simulations are far from definitive: there are countless of ways in which to choose the true ratings and the graph topology.
Our belief that these experiments are indicative of a more general behaviour is due to the peculiarities of the Elo rating systems.
In particular, the Elo update offers diminishing returns: if the rating of $i$ is relatively much larger than the rating of $j$,
if $i$ wins against $j$, the rating of $i$ won't increase by much. This is because the sigmoid function $\sigma$ approaches zero very fast.

Despite this, proving that the maximum rating cannot increase significantly appears much harder.
This is due to the fact that  the maximum rating is not a supermartingale:
if all of the neighbours of the node $i$ with largest true rating have abnormally large rating, the rating of $i$ is likely to increase,
no matter how large it already was.
However, because the zero-sum property of the ratings,
the neighbours of $i$ cannot maintain an abnormally large rating for too long.
Indeed, it appears very unlikely that such balanced but “abnormal” configurations are ever reached.

In the two-player case, the situation is simpler: if player $1$ has Elo rating $x_1$, then player $2$ has $x_2 = -x_1$, by the zero-sum nature.
It is straightforward to check that the Elo ratings are biased towards the true skills:
that is, if player $1$ has rating $x_1$ before a game, then their rating $x'_1$ after the game has $\ex{\abs{x'_1 - \rho_1}} < \abs{x_1 - \rho_1}$, if $\eta < \tfrac12$.

Suppose that the same argument holds for $n > 2$ players:
i.e., typically, Elo ratings are biased towards the true skills.
Biased random walks on $\mbr$ have exponential tails.
So, if $\eta \ll 1/\log n$, then it is super-polynomially unlikely that a given Elo rating will be more than $1$ away from the corresponding true skill, in equilibrium.
A union bound over the $n$ players allows us to deduce that the maximum Elo rating is at most $1$ more than the maximum true skill, with probability at least $1 - 1/n^{10}$.

Unfortunately, we weren't able to make this argument formal for $n > 2$ players.
We leave
proving that indeed, with high probability, the maximum rating remains small for a large number of steps
as an open problem.

§ CONVERGENCE PROOFS: PRELIMINARIES

§ CONVERGENCE PROOFS: BIAS AND VARIANCE

§ CONVERGENCE PROOFS: CURVATURE AND CONCENTRATION

§ CONVERGENCE PROOFS: PARALLEL MATCHES