# FIFA ranking:
Evaluation and path forward

Leszek Szczecinski and Iris-Ioana Roatis L. Szczecinski is with Institut
National de la Recherche Scientifique, Montreal, Canada [e-mail:
<EMAIL_ADDRESS>Roatis is with Imperial College, London, UK
[e-mail<EMAIL_ADDRESS>
###### Abstract

In this work we study the ranking algorithm used by Fédération Internationale
de Football Association (FIFA); we analyze the parameters it currently uses,
show the formal probabilistic model from which it can be derived, and optimize
the latter. In particular, analyzing the games since the introduction of the
algorithm in 2018, we conclude that the game’s “importance” (as defined by
FIFA) used in the algorithm is counterproductive from the point of view of the
predictive capability of the algorithm. We also postulate the algorithm to be
rooted in the formal modelling principle, where the Davidson model proposed in
1970 seems to be an excellent candidate, preserving the form of the algorithm
currently used. The results indicate that the predictive capability of the
algorithm is notably improved by using the home-field advantage (HFA)and the
explicit model for the draws in the game. Moderate, but notable improvement
may be attained by introducing the weighting of the results with the goal
differential, which although not rooted in a formal modelling principle, is
compatible with the current algorithm and can be tuned to the characteristics
of the football competition.

##  1 Introduction

In this work paper we evaluate the algorithm used by Fédération Internationale
de Football Association (FIFA) to rank the international men teams. We also
propose and study simple modifications to improve the prediction capability of
the algorithm.

We are motivated by the fact that the rating and ranking are important
elements of sport competitions and the surrounding entertainment environments.
The rating consists in assigning the team/player a number, often referred to
as “skills” or “strengths”; the ranking is obtained by sorting these numbers
and is also referred to as “power ranking”.

The rating has an informative function providing fans and profane observers
with a quick insight into the relative strength of the teams. For example, the
press is often interested in the “best” teams or the national team reaching
some records position in the ranking.

More importantly, the ranking leads to consequential decisions such as a) the
seeding, i.e., defining which teams play against each other in the
competitions (e.g., used to establish the composition of the groups in the
qualification rounds of the FIFA World Cup), b) the promotion/relegation
(e.g., determining which teams move between the English Premier League (EPL)
and the English Football League Championship, or teams which move between the
Nations Leagues groups), or c) defining the participants in the prestigious
(and lucrative) end-of-the-season competitions (such as Champions League in
European football, Stanley Cup series in National Hockey League (NHL)).

Most of the currently used sport ratings are based on counting of wins/loses
(and draws, when applicable) but in some cases the sport-governing bodies
moved beyond these simple methods and implemented more sophisticated rating
algorithms where the rating levels attributed to the teams are meant to
represent the skills.

In particular, FIFA started a new ranking/rating algorithm in 2018, where the
rating level (skills) assigned to the teams are calculated from the game
outcome, of course, but also from the skills of the teams before the game. The
resulting rating algorithm has a virtue of being simple and defined in a
(mostly) transparent manner.

The main objective of this work is to analyze the FIFA ranking using the
statistical modelling methodology. Considering that the football association
is, by any measure, the most popular sport in the world, it has a value in
itself and follows the line of works which analyzed the past strategies of
ranking used by FIFA, e.g., (Lasek ., 2013), (Ley ., 2019). However, the
approach we propose can be applied to evaluate other algorithms as well, such
as the one used by Fédération Internationale de Volleyball (FIVB), (FIVB,
2020).

In this work we will:

  * •

Derive the FIFA algorithm from the first principles. In particular, we will
define the probabilistic model underlying the algorithm and identify the
estimation method used to estimate the skills.

  * •

Assess the relevance of the parameters used in the current algorithms. In
particular we will evaluate the role played by the change of the adaptation
step according to the game’s importance (as defined by FIFA).

  * •

Optimize the parameters of the proposed model. As a result, we derive an
algorithm which is equally as simple as the FIFA’s one, but allows us to
improve the prediction of the games’ results.

  * •

Propose the modifications of the algorithm which take into account the goals
differential, also known as the margin of victory (MOV). We consider legacy-
compliant algorithms and a new version of the rating.

Our work is organized as follows. In Sec. 2 we describe the FIFA algorithm in
the framework which simplifies the manipulation of the models and the
evaluation of the results. This is also where we clarify the data origin and
make a preliminary evaluation of the relevance of the game’s importance
parameters currently used to control the size of the adaptation step. The
algorithm is then formally derived in Sec. 3 where we also discuss the
evaluation of the results and the batch estimation approach we use. The
incorporation of the MOV into the rating is evaluated in Sec. 4 using two
different strategies. In Sec. 5 we return to the on-line rating, evaluating
and re-optimizing the proposed algorithms, pointing out to the role of the
scale, and commenting on the elements of the FIFA algorithm (the
shootout/knockout rules) which seem to be introduced in an ad-hoc manner for
which the models are not specified. We conclude the work in Sec. 6 summarizing
our findings and in Sec. 6.1 we make en explicit list of recommendations which
may be introduced to improve on the current version of the FIFA algorithm.

##  2 FIFA ranking algorithm

We consider the scenario where there is a total of $M$ teams playing agains
each other in the games indexed with $r=1,\ldots,T$, where $T$ is the number
of games in the observed period. FIFA ranks $M=210$ international teams and,
between June 4, 2018 and  October 16, 2021, there were $T=2964$ games FIFA-
recognized games.

Let us denote the skill of the team $m=1,\ldots,M$ before the game $t$ as
$\theta_{t,m},~{}t\in\mathcal{T}=\\{1,\ldots,T\\}$ which are gathered in a
vector $\boldsymbol{\theta}_{t}=[\theta_{t,1},\ldots,\theta_{t,M}]^{\top}$,
where $(\cdot)^{\top}$ denotes the transpose. The home and the away teams are
denoted by $i_{t}$ and $j_{t}$ respectively.

The game results $y_{t}\in\mathcal{Y}$ are ordinal variables, where the
elements of $\mathcal{Y}=\\{\mathsf{H},\mathsf{D},\mathsf{A}\\}$ represent the
win of the home team ($y_{t}=\mathsf{H}$), the draw ($y_{t}=\mathsf{D}$), and
the win of the away team ($y_{t}=\mathsf{A}$). These ordinal variables are
often transformed into the numerical _scores_
$\check{y}_{t}=\check{y}(y_{t})$: $\check{y}(\mathsf{A})=0$,
$\check{y}(\mathsf{D})=0.5$ and $\check{y}(\mathsf{H})=1$.

The basic rules of FIFA’s rating for a team $m$ which plays in the $t$-th game
are defined as follows

$\displaystyle\theta_{t+1,m}$
$\displaystyle\leftarrow\theta_{t,m}+I_{c_{t}}\delta_{t,m}$ (1)
$\displaystyle\delta_{t,m}$
$\displaystyle=\check{y}_{t,m}-F\big{(}\textstyle\frac{z_{t,m}}{s}\big{)}$ (2)
$\displaystyle F(z)$ $\displaystyle=\frac{1}{1+10^{-z}}$ (3) $\displaystyle
z_{t,m}$ $\displaystyle=\theta_{t,m}-\theta_{t,n}$ (4)

where $s=600$ is the scale111The role of the scale is to ensure that the
values of the skills $\theta_{t,m}$ are situated in a visually comfortable
range; it can be also used when changing the rating algorithm, as is discussed
in Sec. 5.2., $n$ is the index of the team opposing the team $m$ in the $t$-th
game, $\check{y}_{t,m}$ is the “subjective” score of the team $m$ (for the
home team, $m=i_{t}$, $\check{y}_{t,m}=\check{y}_{t}$, and for the away team,
$m=j_{t}$, $\check{y}_{t,m}=1-\check{y}_{t}$). The result produced by the
logistic function, $F(z_{t,m}/s)$ in (2) is referred to as the _expected
score_.

When the team $m$ does not play, its skills do note change, i.e.,
$\theta_{t+1,m}\leftarrow\theta_{t,m}$.

Since $|\delta_{t,m}|\leq 1$, $I_{c_{t}}$ is the maximum allowed update step,
where $c_{t}$ is the game category (or game “importance”) and we can decompose
$I_{c}$ into two components

$\displaystyle I_{c}=K\xi_{c},$ (5)

where $K=5$ and $\xi_{c}$ is a category-dependent adjustment as defined in
Table 1.

$c$ | $I_{c}$ | $\xi_{c}$ | Description | Number
---|---|---|---|---
0 | 5 | 1 | Friendlies outside International Match Calendar windows | 436
1 | 10 | 2 | Friendlies during International Match Calendar windows | 583
2 | 15 | 3 | Group phase of Nations League competitions | 347
3 | 25 | 5 | Play-offs and finals of Nations League competitions | 84
4 | 25 | 5 | Qualifications for Confederations/World Cup finals | 1189
5 | 35 | 7 | Confederation finals up until the QF stage | 209
6 | 40 | 8 | Confederation finals from the QF stage onwards | 52
7 | 50 | 10 | World Cup finals up until QF stage | 56
8 | 60 | 12 | World Cup finals from QF stage onwards | 8
Table 1: Categories, $c$ of the game and the corresponding update steps
$I_{c}=K\xi_{c}$, (FIFA, 2018), where $K=5$ and $\xi_{c}=I_{c}/I_{0}$. The
number of the observed categories between June 4, 2018 and  October 16, 2021
is also given (total number of games is $T=2964$).

The basic equation governing the change of the skills in (2) is next
supplemented with the following rules:

  * •

_Knockout rule_ : in the knockout stage of any competition (which follows the
group stage), instead of (2) we use

$\displaystyle\delta_{t,m}$ $\displaystyle\leftarrow\max\\{0,\delta_{t,m}\\}$
(6)

which guarantees that no points are lost by teams moving out of the group
stage.

  * •

_Shootout rule_ : If the team $m$ wins the game in the shootouts we use

$\displaystyle\check{y}_{t,m}$ $\displaystyle\leftarrow
0.75,\quad\check{y}_{t,n}\leftarrow 0.5,$ (7)

where $n$ is the index of the team which lost.

This rule, however, does not apply in two-legged qualification games if the
shootout is required to break the tie.

We will discuss the effect of the knockout rule later and here we only point
out to the fact that by applying the shootout/knockout we always increase
$\delta_{t,m}$. Thus, while the basic rules (1)-(2) guarantee that the teams
“exchange” the points so the total number of points stays constant, i.e.,
$\sum_{m=1}^{M}\theta_{t,m}=\sum_{m=1}^{M}\theta_{t+1,m}$ (this is a well-
known property of the Elo rating algorithm, (Elo, 1978)), the
shootout/knockout rules increase the total number of points which causes an
“inflation” of the rating. In fact, in the considered period, there were  124
games where the shootout rule, the knockout rule, or both were applied and
this increased the total score by $1739$ points (with the initial total being
$254680$).

The rating we described is published by FIFA since August 2018, roughly on a
per-month basis. The algorithm was initialized on the June 4, 2018, with the
initialization values $\boldsymbol{\theta}_{0}$ based on the previous rating
system.

To run the algorithm, we need to know the initialization
$\boldsymbol{\theta}_{0}$, the presence of conditions which trigger the use of
the knockout/shootout rules, and most importantly, the category/importance of
each game, $c_{t}$. These elements are not officially published so we use here
the unofficial data shown in Football Rankings (2021) which keeps track of the
FIFA rating since June 2018. Using it, we were able to reproduce the ratings
$\boldsymbol{\theta}_{t}$ with a precision of fractions of rating points which
gives us confidence that the categories of the games are assigned according to
the FIFA rules.222Information provided by Football Rankings (2021) is highly
valuable because it is far from straightforward to verify which games are
included in the rating and what their importance $I_{c}$ is. In particular,
the games in the same tournament can be included or excluded from the rating
and in some cases the changes may be done retroactively complicating further
the understanding of the rating results. For example, we had to deal with two
minor exceptions:  • We recognized the victory of Guyana (GUY) over Barbados
(BRB) in the game played on Sept. 6, 2019 already on the date of the game,
while in the FIFA rating, the draw was originally registered and the GUY’s
victory was recognized only later, when BRB was disqualified for having
fielded an ineligible player. • We removed the game Côte d’Ivoire (CIV) vs.
Zambia (ZAM) played on June 19, 2019, where CIV, the winner and ZAM exchanged
2.21 points. The removal of this game from the FIFA-recognized list seems to
be a reason why FIFA changed the ratings of both teams between two official
publications on Dec. 19, 2019 and on Feb. 20, 2020. Namely, the CIV’s rating
was changed from 1380 to 1378 and ZAM’s from 1277 to 1279. This was done
despite both teams not playing at all in this period of time.

Before discussing the suitable models and algorithms we ask a very simple
question: Are the parameters $I_{c}$ defining the “importance” of the game
suitably set? If not, how should we define them to improve the results? The
immediate corollary question is what “improving” the results may mean and, in
general, how to evaluate the results produced by the algorithm.

We note here that the concept of the game importance is not unique to the FIFA
rating and appears also in the FIVB rating, (FIVB, 2020) and the statistical
literature, e.g., (Ley ., 2019, Sec. 2.1.2).

###  2.1 Preliminary evaluation of the FIFA rating

A conventional approach in statistics is to base the performance evaluation on
a metric, called a scoring function, relating the outcome, $y_{t}$ to its
prediction obtained from the estimates at hand (here,
$\boldsymbol{\theta}_{t}$), (Gelman ., 2014).

At this point we want to use only the elements which are clearly defined in
the FIFA ranking and the only explicit predictive element defined in the FIFA
algorithm is the expected score (3),
$F(z_{t}/s)=\mathds{E}[\check{y}_{t}|z_{t}]$, we will base the evaluation on
the metric affected by the mean. Later we will abandon this simplistic
approach.

Using the squared prediction error

$\displaystyle\mathsf{m}(z_{t},y_{t})$
$\displaystyle=\big{(}\check{y}_{t}-F(z_{t}/s)\big{)}^{2}$ (8)

averaged over the large number of games, we obtain the Mean Square Error (MSE)
estimate

$\displaystyle\mathsf{MSE}=\frac{2}{T}\sum_{t=T/2+1}^{T}\mathsf{m}(z_{t},y_{t}),$
(9)

where we use the games in the second half of the observation period to
attenuate the initialization effects. This truncation is somewhat arbitrary of
course but does not affect the results significantly for large $T$.

The MSE in (9) may be treated as an estimate of the expectation,

$\displaystyle\mathsf{MSE}\approx\mathds{E}_{z_{t}}\big{[}\mathds{E}_{y_{t}|z_{t}}[\big{(}\check{y}_{t}-F(z_{t}/s)\big{)}^{2}]\big{]}=\mathds{V}\textrm{ar}[\check{y}_{t}]+\mathds{E}_{z_{t}}\big{[}|B(z_{t},y_{t})|^{2}\big{]},$
(10)

which highlights the bias-variance decomposition, (Duda ., 2001, Ch. 9.3.2)
and where
$\mathds{V}\textrm{ar}[\check{y}_{t}]=\mathds{E}_{z_{t}}\big{[}\mathds{V}\textrm{ar}[\check{y}_{t}|z_{t}]\big{]}$
is the average conditional variance of $\check{y}_{t}$, and
$B(z_{t},y_{t})=F(z_{t}/s)-\mathds{E}[\check{y}_{t}|z_{t}]$ is the estimation
bias of the mean.

Therefore, by reducing the (absolute value of the) bias $B(z_{t},y_{t})$, that
is, by improving the calculation of the expected score $F(z_{t}/s)$, should
manifest itself in a lower value of the MSE, which is calculated as in (9).

Using the MSE, we are now able to assess how the values of the importance
parameters $I_{c}$ (or alternatively, $K$ and $\xi_{c}$) affect the expected
value of the score.

We find the coefficients $K$ and/or $\xi_{c}$ by minimizing the MSE (9) and it
turns out that a simple alternate optimization (one variable $K$ or $\xi_{c}$
is optimized at a time, till convergence) leads efficiently to satisfactory
solutions.333This was done by a line search as we preferred avoiding more
formal, e.g., gradient-based, methods which are not well suited to deal with
the complicated functional relationship resulting from the recursive rating
algorithm.  The results are shown in Table 2 and we observe the following:

  * •

The common update step $K$ increases ten-fold in the optimized setup and it
seems that it is the most important contributor to the improvement of the MSE
(which changes from $\mathsf{MSE}=0.1295$ in the original algorithm to
$\mathsf{MSE}=0.1262$ in the algorithm with fixed-importance games but larger
common adaptation step).

  * •

For the games in the categories well represented in the data, i.e.,
$c\in\\{0,1,2,4,5\\}$, the relative importance of the games $\xi_{c}$ does not
seem to be critically different and for sure does not fall in line with the
values used in the FIFA algorithm. Overall, the optimized $\xi_{c}$ yield a
very small improvement in the MSE comparing to the use of fixed $\xi_{c}$.

In fact that the Friendlies played in the International Match Calendar window
are weighted down ($\xi_{1}=0.6$) comparing to the Friendlies played outside
the window, which is the trend contrary to what the FIFA algorithm does.

  * •

Estimates of $\xi_{c}$ for the categories $c\in\\{3,6,7,8\\}$ should not be
considered as very reliable because the number of games in each of these
categories is rather small (less than $3\%$ of the total). Moreover, the games
in the categories $c=7$ and $c=8$ were observed only in June 2018, during the
2018 World Cup; thus, their effect is most likely very weak in the games from
the second half of the observed batch, see (9), which starts around October
2019.

$\mathsf{MSE}_{\textrm{opt}}$ | $K$ | $\xi_{0}$ | $\xi_{1}$ | $\xi_{2}$ | $\xi_{3}$ | $\xi_{4}$ | $\xi_{5}$ | $\xi_{6}$ | $\xi_{7}$ | $\xi_{8}$
---|---|---|---|---|---|---|---|---|---|---
$0.1295$ | $5$ | $1$ | $2$ | $3$ | $5$ | $5$ | $7$ | $8$ | $10$ | $12$
$0.1262$ | $12$ | $1$ | $2$ | $3$ | $5$ | $5$ | $7$ | $8$ | $10$ | $12$
$0.1262$ | $55$ | $1$ | $1$ | $1$ | $1$ | $1$ | $1$ | $1$ | $1$ | $1$
$0.1250$ | $50$ | $1$ | $0.6$ | $1.8$ | $0.8$ | $1.2$ | $1.1$ | $2.4$ | $0.1$ | $9.9$
Table 2: Parameters $K$ and $\xi_{c}$, in (5), are either fixed (shadowed
cells), or obtained by minimizing the MSE (9). The first row corresponds to
the original FIFA algorithm: $K$ and $\xi_{c}$ are taken from Table 1.

Using a very simple criterion of the MSE derived from the definitions used by
the FIFA algorithm, we obtain results which cast a doubt on the
optimality/utility of the games’ importance parameters, $I_{c}$ proposed by
FIFA.

However, drawing conclusions at this point may be premature. For example,
regarding $K$ (which, after optimization should be much larger than $5$), it
is possible that the relatively short period of observation time ($29$ months)
is not sufficient for small $K$ to guarantee the sufficient convergence but
may pay off in a long run, when smaller values of $K$ will improve the
performance after the convergence is reached. We cannot elucidate this issue
with the data at hand.

On the other hand, to address the concerns regarding the _relative_ importance
weights $\xi_{c}$ the situation is rather different. Even after the
convergence, the weights associated with different game categories should
affect meanigfuly the results. To elucidate this point we will now take a more
formal approach and go back to the “drawing board” to derive the rating
algorithm from the first principles.

##  3 Derivation of the algorithm and batch-rating

To understand and eventually modify the rating algorithm used by FIFA we
propose to cast it in the well defined probabilistic framework. To this end we
define explicitly a model relating the game outcome $y_{t}$ to the skills of
the home-team ($\theta_{i_{t}}$) and the away-team ($\theta_{j_{t}}$), where
the most common assumption is that the probability that, at time $t$, a random
variable $y$ takes the value $y_{t}$, depends on the skills’ difference
$z_{t}=\theta_{t,i_{t}}-\theta_{t,j_{t}}$, i.e.,

$\displaystyle\Pr\left\\{y=y_{t}|\boldsymbol{\theta}_{t}\right\\}$
$\displaystyle=L(z_{t}/s;y_{t})$ (11) $\displaystyle z_{t}$
$\displaystyle=\boldsymbol{x}^{\top}_{t}\boldsymbol{\theta}_{t},$ (12)

where $L(z_{t}/s;y_{t})$ is the _likelihood_ of $\boldsymbol{\theta}_{t}$ (for
a given outcome $y_{t}$) and we define a _scheduling_ vector
$\boldsymbol{x}_{t}=[x_{t,0},\ldots,x_{t,N-1}]^{\top}$ for the game $t$, as

$\displaystyle
x_{t,m}=\mathds{I}\big{[}{i_{t}=m}\big{]}-\mathds{I}\big{[}{j_{t}=m}\big{]},$
(13)

with $\mathds{I}\big{[}{a}\big{]}=1$ when $a$ is true, and
$\mathds{I}\big{[}{a}\big{]}=0$, otherwise. Thus, $x_{t,m}=1$ if the team $m$
is playing at home, $x_{t,m}=-1$ if the team $m$ is visiting, and $x_{t,m}=0$
for all teams $m$ which do not play. This notation allows us to a) deal in a
compact manner with all the skills $\boldsymbol{\theta}_{t}$ for each $t$, and
b) consider the home-field advantage (HFA) or a lack thereof. As before, $s$
is the scale.

We are interested in the on-line rating algorithms, in which the skills of the
participating teams are changed immediately after the results of the game are
known. Nevertheless, we will start the analysis with a batch processing, i.e.,
assuming that the skills $\boldsymbol{\theta}_{t}$ do not vary in time,
$\boldsymbol{\theta}_{t}=\boldsymbol{\theta}$. This is a reasonable approach
if the time window defined $T$ is not too large, so that the skills of the
teams may, indeed, be considered approximately constant. The on-line rating
algorithms will be then derived as approximate solutions to the batch
optimization problem. The purpose of such approach is to a) tie the algorithm
used by FIFA with the theoretical assumptions, which are not spelled out when
the algorithm is presented, b) remove the dependence on the initialization
and/or on the scale, and c) treat the past and present data in the same
manner, e.g., avoiding the partial elimination in the performance metrics, see
(9).

Assuming that the observations are independent when conditioned on the skills,
the rating may be based on the _weighted_ maximum likelihood (ML) estimation
principle

$\displaystyle\hat{\boldsymbol{\theta}}$
$\displaystyle=\mathop{\mathrm{argmin}}_{\boldsymbol{\theta}}\sum_{t\in\mathcal{T}}\xi_{c_{t}}\ell(z_{t}/s;y_{t}),$
(14)

where

$\displaystyle\ell(z_{t}/s;y_{t})$ $\displaystyle=-\log L(z_{t}/s;y_{t}),$
(15)

is a (negated)444The negation in (15) allows us to use a minimization in (14)
which is a very common formulation log-likelihood. The weighting with
$\xi_{c_{t}}\in(0,1]$ is used in the estimation literature to take care of the
outcomes which are more or less reliable, (Hu Zidek, 2001), (Amiguet, 2010).
In our problem, the reliability is associated with the game category, $c_{t}$,
so $\xi_{c}$ denotes the weight of the category $c$. Since multiplication of
$\xi_{c}$ by a common factor is irrelevant for minimization, we fix
$\xi_{0}=1$.

We may solve (14) using the steepest descent

$\displaystyle\hat{\boldsymbol{\theta}}\leftarrow\hat{\boldsymbol{\theta}}-\mu/s\sum_{t}\boldsymbol{x}_{t}\xi_{c_{t}}g(z_{t}/s;y_{t}),$
(16)

where $\mu$ is the adaptation step and

$\displaystyle g(z;y)$
$\displaystyle=\frac{\,\mathrm{d}}{\,\mathrm{d}z}\ell(z;y).$ (17)

The on-line version of (16) is obtained replacing the batch-optimization with
the stochastic gradient (SG) which updates the solution each time a new
observation becomes available, i.e.,

$\displaystyle\boldsymbol{\theta}_{t+1}\leftarrow\boldsymbol{\theta}_{t}-K\xi_{c_{t}}\boldsymbol{x}_{t}g(z_{t}/s;y_{t}),$
(18)

where the update amplitude is controlled by the weight $\xi_{x_{t}}$ and the
step $K$ which absorbs the scale $s$.

###  3.1 Davidson model and Elo algorithm

The rating depends now on the choice of the likelihood function $L(z;y)$ and
we opt here for the Davidson model, (Davidson, 1970), being a particular case
of the multinomial model used also in Egidi Torelli (2021)

$\displaystyle L(z;\mathsf{H})$ $\displaystyle=\frac{10^{0.5(z+\eta
b)}}{10^{0.5(z+\eta b)}+\kappa+10^{-0.5(z+\eta b)}},$ (19) $\displaystyle
L(z;\mathsf{A})$ $\displaystyle=P(-z;\mathsf{H}),$ (20) $\displaystyle
L(z;\mathsf{D})$ $\displaystyle=\kappa\sqrt{L(z;\mathsf{H})L(z;\mathsf{A})},$
(21)

where $\eta$ is a HFA modelling the apparent increase in the skills of the
local team, the indicator $b=\mathds{I}\big{[}{\textrm{game is played in the
home-team country}}\big{]}$ allows us to distinguish between the games played
on the home or the neutral venues,555Out of $T=2964$ games we considered,
$768$ were played on neutral venues. To verify the venue we used (The Roon Ba,
2021) and (SoccerWay, 2021). and $\kappa$ controls for the presence of the
draws.

The choice of this model is motivated by the fact that it leads to a simple
algorithmic update of the skills generalizing the Elo rating algorithm,
(Szczecinski Djebbi, 2020). It also becomes equivalent to the latter for a
particular value of $\kappa=2$ and $\kappa=0$. These relationships make
possible (or at the least – ease) the comparison with the rating algorithm
currently used by FIFA which is also based on the Elo algorithm.

Using (19)-(21) in (17), with straightforward algebra we obtain

$\displaystyle g(z;y)$
$\displaystyle=\frac{\,\mathrm{d}}{\,\mathrm{d}z}\ell(z;y)$ (22)
$\displaystyle=-\ln 10(\check{y}-F_{\kappa}(z)),$ (23)

where $\check{y}$ is the “score” of the game which we already defined, and

$\displaystyle F_{\kappa}(z)$
$\displaystyle=\frac{\frac{1}{2}\kappa+10^{0.5(z+\eta b)}}{10^{0.5(z+\eta
b)}+\kappa+10^{-0.5(z+\eta b)}}$ (24)

has the meaning of the conditional expected score,
$F_{\kappa}(z)=\mathds{E}[\check{y}|z]=\sum_{y\in\mathcal{Y}}\check{y}L(z;y)$.

Therefore, the SG algorithm (18) becomes

$\displaystyle\boldsymbol{\theta}_{t+1}\leftarrow\boldsymbol{\theta}_{t}+K\xi_{t}\boldsymbol{x}_{t}\big{(}\check{y}_{t}-F_{\kappa}(z_{t}/s)\big{)}$
(25)

and it obviously has the form of the Elo and FIFA rating algorithms, see
(1)-(3), except that we use $F_{\kappa}(z)$ while the former use $F(z)$. Note
that the step $K$ in (25) absorbs the term $\ln 10$ from (23).

It is easy to see that for $\eta=0$ and $\kappa=0$ (i.e., when
$L(z;\mathsf{D})=0$ and the draws are ignored) we have $F_{0}(z)=F(z)$ which
is simply a logistic function as in the Elo algorithm. Furthermore, for
$\eta=0$ and $\kappa=2$ we obtain $F_{2}(z)=F(z/2)$ and thus (25) is again
equivalent to the Elo rating algorithm but with the doubled scale value.

While we conclude that the FIFA rating algorithm may be seen as the instance
of the maximum weighted likelihood estimation, this is, of course, a “reverse-
engineered” hypothesis because the FIFA document, (FIFA, 2018), does not
mention any remotely similar concept.

###  3.2 Regularized batch rating

In order to go beyond the limitation of the SG optimization and to avoid the
problems related to the removal of the significant portion of the data (meant
to eliminate the initialization effects during evaluation, see (9)) we may
focus on the original problem defined in (14) for the entire set of data. That
is, we ignore now the on-line rating aspect and rather focus on the evaluation
of the model and the optimization criterion that underlie the algorithm.

We start noting that the problem (14) is, in general, ill-posed: since the
solution depends only on the differences between the skills, $z_{t}$, all
solutions $\hat{\boldsymbol{\theta}}$ and
$\hat{\boldsymbol{\theta}}+\theta_{\textrm{o}}\boldsymbol{1}$ are equivalent
because the differences $z_{t}$ are independent from “origin” value
$\theta_{\textrm{o}}$. To remove this ambiguity we may _regularize_ the
problem as

$\displaystyle\hat{\boldsymbol{\theta}}$
$\displaystyle=\mathop{\mathrm{argmin}}_{\boldsymbol{\theta}}J(\boldsymbol{\theta})$
(26) $\displaystyle J(\boldsymbol{\theta})$
$\displaystyle=\sum_{t\in\mathcal{T}}\xi_{c_{t}}\ell(z_{t}/s;y_{t})+\frac{\alpha}{2s^{2}}\|\boldsymbol{\theta}\|^{2},$
(27)

where $\alpha$ is the regularization parameter and we have opted for a so-
called ridge regularization (Hastie ., 2009, Ch. 3.4.1).

Under the model (19)-(21), the regularized batch-optimization problem (26) is
useful to resolve another difficulty. Namely, if there is a team $m$ having
registered only wins, i.e., when $\forall i_{t}=m,~{}y_{t}=\mathsf{H}$ and
$\forall j_{t}=m,~{}y_{t}=\mathsf{A}$, then (14) cannot be solved (or rather,
$\hat{\theta}_{m}\rightarrow\infty$) because $J(\boldsymbol{\theta})$ does not
limit the value of $\theta_{m}$. Such a solution not only is unattainable
numerically but is, in fact meaningless and the regularization (26) settles
this issue.666The same problem arises, of course, when a team registers a
sequences of pure losses. This is not a hypothetical issue and in the official
FIFA games three teams registered the streaks of unique wins or losses
(without any other results): Tonga (three wins), Eritrea (two losses), and
American Samoa (four losses). Thus, the attempt to solve the batch-
optimization problem without regularization (i.e., with $\alpha=0$) would
yield $\hat{\theta}_{m}=\infty$, for $m$ being index of Tonga.

The estimated skills, $\hat{\boldsymbol{\theta}}$ depend now on the weights
$\xi_{c}$, on the regularization parameter $\alpha$, and on the model
parameters $\eta$ and $\kappa$. If unknown, all these parameters must be
optimized.

As for the optimization criterion, we recall that the FIFA algorithm only
specified the expected score, so the quadratic error (8) was allowed us to
evaluate the algorithm and stay within the boundaries of its definitions. Now,
however, with the explicit skills-outcome model, we may go beyond this
limitation and may use the prediction metrics known in the machine learning
such as the (negated) log-score, (Gelman ., 2014)

$\displaystyle\mathsf{m}^{\textrm{ls}}(z_{t};y_{t})$
$\displaystyle=\ell(z_{t}/s;y_{t}),$ (28)

often preferred due to its compatibility with the log-likelihood used as the
optimization criterion, or the accuracy score, (Lasek Gagolewski, 2020)

$\displaystyle\mathsf{m}^{\textrm{acc}}(z_{t};y_{t})$
$\displaystyle=\mathds{I}\big{[}{y_{t}=\mathop{\mathrm{argmax}}_{y}L(z_{t}/s;y)}\big{]},$
(29)

which equals one if the event with the largest predicted probability was
actually observed, otherwise it is zero.

Furthermore, thanks to the batch-rating we are able to consider the entire
data set in the performance evaluation by averaging the scoring function (28)
or (29) over all games

$\displaystyle\mathsf{LS}$
$\displaystyle=\frac{1}{T}\sum_{t\in\mathcal{T}}\mathsf{m}^{\textrm{ls}}\big{(}\boldsymbol{x}_{t}^{\top}\hat{\boldsymbol{\theta}}_{\backslash{t}},y_{t}\big{)},$
(30) $\displaystyle\mathsf{ACC}$
$\displaystyle=\frac{1}{T}\sum_{t\in\mathcal{T}}\mathsf{m}^{\textrm{acc}}\big{(}\boldsymbol{x}_{t}^{\top}\hat{\boldsymbol{\theta}}_{\backslash{t}},y_{t}\big{)},$
(31)

where

$\displaystyle\hat{\boldsymbol{\theta}}_{\backslash{t}}$
$\displaystyle=\mathop{\mathrm{argmin}}_{\boldsymbol{\theta}}J_{\backslash{t}}(\boldsymbol{\theta}),$
(32) $\displaystyle J_{\backslash{t}}(\boldsymbol{\theta})$
$\displaystyle=\sum_{\begin{subarray}{c}l\in\mathcal{T}\\\ l\neq
t\end{subarray}}\xi_{c_{l}}\ell(\boldsymbol{x}_{l}^{\top}\boldsymbol{\theta}/s;y_{l})+\frac{\alpha}{2s^{2}}\|\boldsymbol{\theta}\|^{2}.$
(33)

In plain words, for given parameters ($\alpha$, $\kappa$, $\eta$, $\xi_{c}$),
we find the skills $\hat{\boldsymbol{\theta}}_{\backslash{t}}$ from all, but
the $t$-th game [this is (32)-(33)], and next use them to predict the results
$y_{t}$; we repeat it for all $t\in\mathcal{T}$, summing the obtained scores.
This is the well-known leave-one-out (LOO) cross-validation strategy (Hastie
., 2009, Sec. 2.9), (Duda ., 2001, Ch. 9.6.2): no data is discarded when
calculating the metrics (30)-(31) and this comes with the price of having to
find $\hat{\boldsymbol{\theta}}_{\backslash{t}}$ for all $t\in\mathcal{T}$. To
diminish the computational load, we opt here for the approximate leave-one-out
(ALO) cross-validation (Rad Maleki, 2020) based on the local quadratic
approximation of the optimization function defined for all the data. Details
are given in Appendix A.

Although both, the average log-score in (30) and the accuracy (31) can be now
optimized with respect to $\alpha$, $\kappa$, $\eta$, and/or $\xi_{c}$, we
only optimize the log-score whose optimal value is denoted as
$\mathsf{LS}_{\textrm{opt}}$; the resulting accuracy, $\mathsf{ACC}$ will be
also shown. It is, of course, possible to optimize the log-score with respect
to any subset of parameters.

Again, we used the alternated minimization: $\mathsf{LS}$ was minimized with
respect to one parameter at a time: $\alpha$, $\kappa$, $\eta$, or $\xi_{c}$,
till no improvement was observed. This simple strategy led to the minimum
$\mathsf{LS}_{\textrm{opt}}$ which turned out to be independent of various
starting points we used.777Although we cannot prove the solution to be global,
in all our observations the log-score functions seemed to be unimodal.

A quick comment may be useful regarding the interpretation of the performance
metrics. The accuracy (31) is easily understandable: it is an average number
of the events which were predicted correctly (as those with the maximum
likelihood $L(z_{t}/s;y)$). On the other hand, the metric (30) may be
represented as $\exp(-\mathsf{LS})=[\prod_{t=1}^{T}L(z_{t}/s;y_{t})]^{1/T}$
which is a geometric mean of the predicted probabilities assigned to the
events which were actually observed. While the accuracy metric penalizes the
wrong guesses with zero (so $\mathsf{ACC}\in[0,1]$), the log-score penalizes
them via the logarithmic function, which may be arbitrarily large (so
$\mathsf{LS}\in(0,\infty)$).

However, the fundamental difference between the two metrics is that we can use
the accuracy without specifying the distribution for all possible outcomes but
we cannot calculate the log-score in such a case.888 The common confusion is
to interpret the function $F(z_{t}/s)$ in the Elo/FIFA algorithm as the
probability of the home win, and the value $1-F(z_{t}/s)$, as the probability
of an away win. This, of course, implies that the draw probability equals
zero. With such an interpretation, we can still calculate the accuracy metric
even if we never predict the draw. On the other hand we cannot calculate the
log-score, because when the draw occurs, we have undefined metric
$\mathsf{m}^{\textrm{ls}}(z_{t}/s;\mathsf{D})\rightarrow\infty$.

$\mathsf{LS}_{\textrm{opt}}$ | $\alpha$ | $\eta$ | $\kappa$ | $\xi_{0}$ | $\xi_{1}$ | $\xi_{2}$ | $\xi_{3}$ | $\xi_{4}$ | $\xi_{5}$ | $\xi_{6}$ | $\xi_{7}$ | $\xi_{8}$ | $\mathsf{ACC}[\%]$
---|---|---|---|---|---|---|---|---|---|---|---|---|---
$0.960$ | $1.7$ | $0$ | $2.0$ | $1$ | $2.0$ | $3.0$ | $5.0$ | $5.0$ | $7.0$ | $8.0$ | $10.0$ | $12.0$ | $55$
$0.948$ | $0.2$ | $0$ | $2.0$ | $1$ | $1.0$ | $1.0$ | $1.0$ | $1.0$ | $1.0$ | $1.0$ | $1.0$ | $1.0$ | $56$
$0.948$ | $0.3$ | $0$ | $2.0$ | $1$ | $0.9$ | $0.7$ | $0.8$ | $0.9$ | $1.2$ | $0.7$ | $0.8$ | $1.1$ | $55$
$0.918$ | $0.3$ | $0.4$ | $2.0$ | $1$ | $1.0$ | $1.0$ | $1.0$ | $1.0$ | $1.0$ | $1.0$ | $1.0$ | $1.0$ | $56$
$0.860$ | $0.4$ | $0.3$ | $0.8$ | $1$ | $1.0$ | $1.0$ | $1.0$ | $1.0$ | $1.0$ | $1.0$ | $1.0$ | $1.0$ | $61$
$0.860$ | $0.5$ | $0.3$ | $0.8$ | $1$ | $0.8$ | $0.7$ | $0.8$ | $1.0$ | $1.2$ | $0.8$ | $1.0$ | $1.0$ | $61$
Table 3: Batch-rating parameters obtained via minimization of the log-score
(30). The parameters ($\alpha$, $\kappa$, $\eta$, $\xi_{c}$) are either fixed
(shadowed cells), or obtained via optimization. The upper-part results
correspond to the conventional FIFA algorithm: using $\kappa=2$ and $\eta=0$,
the expected score is calculated using a logistic function.

The results obtained are shown in Table 3 and indicate that

  * •

The data does not provide evidence for using a category-dependent weights
$\xi_{c}$. There is actually a slight indication that the optimal weights of
the Friendlies within the IMC (category $c=1$) and Group phase of Nations
Leagues (category $c=2$) are slightly _smaller_ than the weight of the regular
Friendlies. This stands in contrast to the FIFA algorithm which doubles the
weight $\xi_{1}$ of the Friendlies played in the IMC and triples the weight of
$\xi_{2}$.

In fact, the results obtained using the FIFA weights $\xi_{c}$ are _worse_
than those obtained using constant weights $\xi_{c}=1$ (i.e., essentially
ignoring the possibility of weighting). Even with the argument of having a
small number of games in some categories (such as a World Cup), it is very
unlikely that observing more games will speak in favor of variable weights and
almost surely not in favor of the highly disproportionate weights used in the
FIFA algorithm.

  * •

A notable improvement in the prediction capacity as measured by the log-score
is obtained by considering the HFA. The value $\eta\in\\{0.3,0.4\\}$ emerges
from the optimization fit and we note that $\eta=0.25$ was used in
eloratings.net (2020).999Therein, the unnormalized value $\eta s=100$ is
reported and since $s=400$, we obtain $\eta=0.25$.

  * •

A more important improvement is obtained by optimizing the parameter $\kappa$
which takes into account the draws and their frequency as discussed in
Szczecinski Djebbi (2020).

It is interesting to compare the parameters found by optimization with the
simplified formulas proposed in Szczecinski Djebbi (2020, Sec. 3.2)

$\displaystyle\eta$
$\displaystyle=\log_{10}\frac{f_{\mathsf{H}}}{f_{\mathsf{A}}}$ (34)
$\displaystyle\kappa$
$\displaystyle=\frac{f_{\mathsf{D}}}{\sqrt{f_{\mathsf{H}}f_{\mathsf{A}}}}\approx\frac{2f_{\mathsf{D}}}{1-f_{\mathsf{D}}}.$
(35)

where $f_{y},y\in\mathcal{Y}$ are empirical frequencies of outcomes. We can
consider separately the games played on the neutral venues and calculate these
frequencies as $f^{\textrm{neut.}}_{\mathsf{A}}=0.37$,
$f^{\textrm{neut.}}_{\mathsf{D}}=0.24$,
$f^{\textrm{neut.}}_{\mathsf{H}}=0.39$, and those played on home venues as
$f^{\textrm{hfa}}_{\mathsf{A}}=0.27$, $f^{\textrm{hfa}}_{\mathsf{D}}=0.22$,
$f^{\textrm{hfa}}_{\mathsf{H}}=0.51$, which yields

$\displaystyle\kappa^{\textrm{hfa}}$ $\displaystyle=0.61$
$\displaystyle\eta^{\textrm{hfa}}$ $\displaystyle=0.28$ (36)
$\displaystyle\kappa^{\textrm{neut.}}$ $\displaystyle=0.63$
$\displaystyle\eta^{\textrm{neut.}}$ $\displaystyle=0.02.$ (37)

The parameter $\eta^{\textrm{hfa}}$ predicted by (34) is practically equal to
the one obtained by optimization. And while the parameters
$\kappa^{\textrm{hfa}}$ and $\kappa^{\textrm{neut.}}$ are slightly different
from the one predicted by (35), using them in the rating, we obtained
$\mathsf{LS}_{\textrm{opt}}=0.868$, which is still notably better than using
the conventional FIFA rating. This is interesting because finding the
parameters $\eta$ and $\kappa$ from the frequencies of the games not only
avoids optimization but also provides a simple empirical justification.

##  4 Margin of victory

In the search for a possible improvement of the rating we want to consider now
the use of the MOV variable, defined by the difference of the goals scored by
each team, and denoted by $d_{t}$. With that regard, the most recent works
adopt two conceptually different approaches.

The first one keeps the structure of the known rating algorithm (such as FIFA
algorithm) and modifies it by changing the adaptation step size as a function
of $d_{t}$. This was already done in eloratings.net (2020), Hvattum Arntzen
(2010), Silver (2014), Ley . (2019), and Kovalchik (2020), and is conceptually
similar to the weighting according to the game-category we consider in the
previous section.

Second approach changes the model relating the skills to the MOV variable
$d_{t}$ and was already studied before in Maher (1982), Ley . (2019), Lasek
Gagolewski (2020), Szczecinski (2020). We will focus on the simple proposition
from Lasek Gagolewski (2020) building on the formulation of Karlis Ntzoufras
(2008).

###  4.1 MOV via weighting

For the context, we show in Table 4 the number of games and their percentage
of the total, depending on the value of the MOV variable $d$. While, in
principle, it is possible to use directly $d$, it is customary to consider
their absolute value, $|d|$.

The Elo/FIFA algorithms (25) can be easily modified as follows, to take the
MOV variable into account:

$\displaystyle K_{c,d}=K\xi_{c}\zeta_{d},$ (38)

where, as before, $K$ is the common step, $\xi_{c}$ is the weight associated
with the game-category $c$, and $\zeta_{d}$ is the function of the MOV-
variable $d$.

$|d|=0$ | $|d|=1$ | $|d|=2$ | $|d|=3$ | $|d|=4$ | $|d|=5$ | $|d|>6$
---|---|---|---|---|---|---
$678$ | $1070$ | $543$ | $337$ | $162$ | $77$ | $97$
$22\%$ | $36\%$ | $18\%$ | $11\%$ | $6\%$ | $3\%$ | $2\%$
Table 4:  Number of games till October 15, 2021 which finished with the goal
difference $|d|$ (the fractions are not adding to 100% due to rounding).

For example, (eloratings.net, 2020) uses

$\displaystyle\zeta_{d}$ $\displaystyle=\begin{cases}1&|d|\leq 1\\\
1.5&|d|=2\\\ 1.75+0.125(|d|-3)&|d|\geq 3\end{cases}.$ (39)

The similar propositions may be found in Hvattum Arntzen (2010) (in the
context of association football), in Kovalchik (2020) (to rate the tennis
players), or in Silver (2014) (for rating of the teams in American football).

To elucidate how useful such heuristics are, we note that the problem is very
similar to the importance-weighting we analyzed before; the difference resides
in the fact that the weighting depends now on the product $\xi_{c}\zeta_{d}$.
We may thus reuse our optimization strategy to find the optimal weights for
the games with different values of $|d|$.

To this end we discretize $|d|$ into $V+1$ MOV-categories, $v=0,\ldots,V$ and
we use a very simple mapping $v=|d|$ for $v<V$ and $v=V\iff|d|\geq V$. For
example, with $V=2$, $\zeta_{0}$ weights the draws ($|d|=0$), $\zeta_{1}$
weights the games with one goal difference ($|d|=1),v=0,1$ and $\zeta_{2}$
weights the games with more than one goal difference ($|d|\geq 1$).

Breaking with the predefined functional relationship as the one shown in (39)
we are more general than the latter, e.g., treating the cases $|d|=0$ and
$|d|=1$ separately. This makes sense since, not only they are the most
frequent events, corresponding, respectively, to $22\%$ and $36\%$ of the
total, see Table 4, but also they correspond to the events of draw and
win/loss treated differently by the algorithm.

On the other hand, we are also less general due to the merging of the events
$|d_{t}|\geq V$, although this effect will decrease with $V$, simply because
there will be very few observations as may be understood from Table 4. For
example, with $V=4$, the weighting $\zeta_{4}$ will be the same for the events
with $|d|=4$ and $|d|>4$ but the latter make only $5\%$ of the total.

We consider again the game categories defined in Table 1 and thus we solve now
the problem

$\displaystyle\hat{\boldsymbol{\theta}}$
$\displaystyle=\mathop{\mathrm{argmin}}_{\boldsymbol{\theta}}\sum_{t\in\mathcal{T}}\xi_{c_{t}}\zeta_{v_{t}}\ell(z_{t}/s;y_{t})+\frac{\alpha}{2s^{2}}\|\boldsymbol{\theta}\|^{2},$
(40)

where $v_{t}$ is the index of the MOV variable $d_{t}$. To remove ambiguity of
the solution, we set $\xi_{0}=1$ and $\zeta_{0}=1$.

The parameters $\xi_{c}$, $\zeta_{v}$, $\eta$, $\kappa$, and $\alpha$ will be
chosen again using the ALO approach we described in Sec. 3.2, that is, by
optimizing the log-score criterion (30). The results shown in Table 5 allow us
to conclude that:

  * •

The optimization of the MOV-weights $\zeta_{v}$ (while keeping $\xi_{c}=1$)
yields $\mathsf{LS}_{\textrm{opt}}=0.937$ and the optimization of $\xi_{c}$
(with $\zeta_{v}=1$) yields $\mathsf{LS}_{\textrm{opt}}=0.948$ (see Table 3).
By comparing them, we see that weighting of the MOV-categories is more
beneficial than weighting of the game-categories. Therefore, there is little
improvement in considering the category-related weights, $\xi_{c}$.

  * •

The optimization indicates that $\zeta_{v}$ defined by (39) is suboptimal. In
particular, the optimal MOV weights, $\zeta_{v}$ are monotonically growing (as
foreseen by the heuristics) only for $|d|\geq 1$ and the draws (i.e.,
$|d_{t}|=0$) have a weight which is more important that the weights of the
events $|d_{t}|=1$; thus, these two events should not be merged together, nor
we should impose a particular functional form for the weights.

  * •

The best improvement in the prediction is obtained again by optimizing the
parameters $\eta$ and $\kappa$ of the Davidson model together with the MOV
weights $\zeta_{v}$.

$\mathsf{LS}_{\textrm{opt}}$ | $V$ | $\alpha$ | $\eta$ | $\kappa$ | $\boldsymbol{\xi}$ | $\zeta_{0}$ | $\zeta_{1}$ | $\zeta_{2}$ | $\zeta_{3}$ | $\zeta_{4}$ | $\zeta_{5}$ | $\zeta_{6}$ | $\mathsf{ACC}[\%]$
---|---|---|---|---|---|---|---|---|---|---|---|---|---
$0.949$ | $6$ | $0.9$ | $0$ | $2.0$ | $\boldsymbol{1}$ | $1$ | $1$ | $1.5$ | $1.75$ | $1.875$ | $2.0$ | $4.2$ | $56$
$0.937$ | $6$ | $0.2$ | $0$ | $2.0$ | $\boldsymbol{1}$ | $1$ | $0.3$ | $0.5$ | $0.7$ | $1.0$ | $1.3$ | $2.4$ | $55$
$0.935$ | $6$ | $0.2$ | $0$ | $2.0$ | $\hat{\boldsymbol{\xi}}$ | $1$ | $0.3$ | $0.5$ | $0.7$ | $1.0$ | $1.5$ | $2.3$ | $55$
$0.906$ | $6$ | $0.2$ | $0.4$ | $2.0$ | $\boldsymbol{1}$ | $1$ | $0.3$ | $0.5$ | $0.7$ | $1.0$ | $1.6$ | $3.0$ | $56$
$0.852$ | $6$ | $0.3$ | $0.3$ | $0.8$ | $\boldsymbol{1}$ | $1$ | $0.3$ | $0.5$ | $0.7$ | $1.0$ | $1.6$ | $3.0$ | $62$
$0.853$ | $4$ | $0.3$ | $0.3$ | $0.8$ | $\boldsymbol{1}$ | $1$ | $0.3$ | $0.5$ | $0.8$ | $1.4$ | $\times$ | $\times$ | $62$
$0.854$ | $2$ | $0.2$ | $0.3$ | $0.8$ | $\boldsymbol{1}$ | $1$ | $0.3$ | $0.8$ | $\times$ | $\times$ | $\times$ | $\times$ | $62$
$0.857$ | $1$ | $0.2$ | $0.3$ | $0.8$ | $\boldsymbol{1}$ | $1$ | $0.6$ | $\times$ | $\times$ | $\times$ | $\times$ | $\times$ | $61$
Table 5: Batch-rating parameters obtained via minimization of the log-score
(30) with weighting of the MOV-variables. The parameters ($\alpha$, $\kappa$,
$\eta$, $\boldsymbol{\xi}$, $\zeta_{v}$) are either fixed (shadowed cells), or
obtained through optimization; to save space, in the sole case when the
parameters $\xi_{c}$ are optimized, their optimal values are gathered in the
vector $\hat{\boldsymbol{\xi}}=[1.0,1.0,0.9,1.2,1.1,1.8,1.9,1.4,5.6]$.

###  4.2 MOV via modelling

A different approach to deal with the MOV relies on the integration of the
latter into the formal model relating the skills $\boldsymbol{\theta}_{t}$ and
the observed MOV variable $d_{t}$.

A simple approach proposed in Karlis Ntzoufras (2008) relies on a direct
modelling of the goal difference using the Skellam’s distribution

$\displaystyle\Pr\left\\{d_{t}=d|\boldsymbol{\theta}_{t}\right\\}$
$\displaystyle=L(z_{t};d_{t})$ (41)
$\displaystyle=\mathrm{e}^{-(\mu_{\textrm{h},t}+\mu_{\textrm{a},t})}\left(\frac{\mu_{\textrm{h},t}}{\mu_{\textrm{a},t}}\right)^{d/2}I_{|d|}(2\sqrt{\mu_{\textrm{h},t}\mu_{\textrm{a},t}}),$
(42)

where $I_{v}(t)$ is the modified Bessel function of order $v$ and
$\mu_{\textrm{h},t}$ and $\mu_{\textrm{h},t}$ are means of the Poisson
variables modelling the home- and away- goals. The latter are functions of the
skills’ difference $z_{t}$, (Karlis Ntzoufras, 2008, Sec. 2.2)

$\displaystyle\mu_{\textrm{h},t}$
$\displaystyle=\mathrm{e}^{c+z_{t}+b\eta},\quad\mu_{\textrm{a},t}=\mathrm{e}^{c-z_{t}-b\eta},$
(43)

where is $c$ is a constant and, as before, $\eta$ is the HFA
coefficient.101010For the home team we add– and for the away team – subtract
$b\eta$ in the exponent. This is different from (Karlis Ntzoufras, 2008),
(Lasek Gagolewski, 2020), where only the home team benefits from the HFA boost
while the away team is not penalized, see (Karlis Ntzoufras, 2008, Eq.
(2.2)-(2.3)). Of course, we can rewrite (43) as
$\mu_{\textrm{h},t}=\mathrm{e}^{c^{\prime}+z_{t}+b\eta^{\prime}}$,
$\mu_{\textrm{a},t}=\mathrm{e}^{c^{\prime}+z_{t}}$ with $c^{\prime}=c-b\eta$
and $\eta^{\prime}=2\eta$ but it makes sense only when the HFA is always
present, as then $b\eta=\eta$. While this condition holds in the context of
football leagues considered in Karlis Ntzoufras (2008) and in Lasek Gagolewski
(2020), this is not the case in the international FIFA games, which can be
played on the neutral venues.

The model (41) is a particular case of a more general form shown in Karlis
Ntzoufras (2008), which allowed us to model the offensive and the defensive
skills. Here, however, we are interested in rating and thus one skill per team
should be used. As noted in Ley . (2019), Lasek Gagolewski (2020) this offers
a sufficient prediction capability avoiding the problem of over-
parametrization due to doubling of the number of skills.

Using (43) in (41), the following log-likelihood is obtained

$\displaystyle\ell(z;d)$ $\displaystyle=-\log L(z;d)$ (44)
$\displaystyle=(\mu_{\textrm{h}}+\mu_{\textrm{a}})-d(z+b\eta)-2\mathrm{e}^{c}-\log\tilde{I}_{|d|}(2\mathrm{e}^{c})$
(45)

where, for numerical stability it is convenient to use an exponentially
modified form of the Bessel function,
$\tilde{I}_{v}(t)=I_{v}(t)\mathrm{e}^{-t}$, available in many computation
packages.

The derivative of (44) is given by

$\displaystyle g(z;d)$
$\displaystyle=\frac{\,\mathrm{d}}{\,\mathrm{d}z}\ell(z;d)=-(d-\overline{F}(z)),$
(46) $\displaystyle\overline{F}(z)$
$\displaystyle=\mu_{\textrm{h}}-\mu_{\textrm{a}}=\mathrm{e}^{c}(\mathrm{e}^{z+b\eta}-\mathrm{e}^{-z-b\eta}).$
(47)

The batch rating consists then in solving the following problem:

$\displaystyle\hat{\boldsymbol{\theta}}$
$\displaystyle=\mathop{\mathrm{argmin}}_{\boldsymbol{\theta}}\sum_{t\in\mathcal{T}}\ell(z_{t}/s;d_{t})+\frac{\alpha}{2s^{2}}\|\boldsymbol{\theta}\|^{2}$
(48)

and the SG implementation of the ML principle will produce the algorithm

$\displaystyle\boldsymbol{\theta}_{t+1}\leftarrow\boldsymbol{\theta}_{t}+K\boldsymbol{x}_{t}\big{(}d_{t}-\overline{F}(z_{t}/s)\big{)},$
(49)

which is again written in a form similar to the FIFA rating algorithm, where
the goal difference $d_{t}$ plays the role of the “score”, and
$\overline{F}(z_{t}/s)=\mathds{E}[d_{t}|z_{t}]$ is the expected score. The
algorithm (49) can be also obtained by applying the Poisson model to the goals
scored by each of the teams (Lasek Gagolewski, 2020).

To calculate the log-score, we have to merge the events $d<0$ (away-win) and
$d>0$ (home-win). Since the closed-form formulas do not exist we do it
approximately by truncated sums

$\displaystyle\mathsf{m}^{\textrm{ls}}(z;\mathsf{A})$
$\displaystyle=-\log\sum_{d=-D}^{-1}L(z;d),\quad\mathsf{m}^{\textrm{ls}}(z;\mathsf{D})=-\log
L(z;0),\quad\mathsf{m}^{\textrm{ls}}(z;\mathsf{H})=-\log\sum_{d=1}^{D}L(z;d),$
(50)

where we used $D=50$ which guaranteed that
$|1-\sum_{d=-D}^{D}L(z;d)|<10^{-4}$.

The results shown in Table 6 indicate that, with this very simple approach
(with only two parameters of the model which must be optimized) we are able to
improve over the MOV-weighting strategy and this should be attributed to the
use of a formal skills-outcome model. The price to pay for the improvement
lies in the change of the entire algorithm and in abandoning of the legacy of
the Elo algorithm.

Moreover, the possible implementation issues may arise since the expected
score (47) is theoretically unbounded. Whether the improvement of the log-
score from $\mathsf{LS}=0.857$ (in the MOV-weighting, see Table 5) to
$\mathsf{LS}=0.845$ in the Skellam’s MOV model is worth the change and the
implementation risks, is at least debatable.

$\mathsf{LS}_{\textrm{opt}}$ | $\alpha$ | $\eta$ | $c$ | $\mathsf{ACC}[\%]$
---|---|---|---|---
$0.845$ | $0.21$ | $0.20$ | $0$ | $61$
Table 6: Batch-rating parameters obtained via minimization of the log-score
(30) using the Skellam’s model (44).

##  5 On line rating

original | no shootouts | no knockouts | no shootouts
---|---|---|---
| rules | rules | /knockouts rules
BEL (1832.3) | BEL (1831.0) | BRA (1775.9) | FRA (1768.4)
BRA (1820.4) | BRA (1817.5) | FRA (1770.1) | BRA (1767.7)
FRA (1779.2) | FRA (1778.2) | BEL (1759.2) | BEL (1757.4)
ITA (1750.5) | ENG (1740.2) | ITA (1730.5) | ITA (1711.1)
ENG (1750.2) | ITA (1733.0) | ENG (1711.9) | ENG (1701.3)
Table 7: Ranking of the top teams: Belgium (BEL), Brazil (BRA), England (ENG),
France (FRA), and Italy (ITA). The original FIFA algorithm and its modified
rules are considered.

Before starting a metrics-based comparison of the on-line algorithms, in Sec.
5.1 we will address the use of the knockout/shootout rules (6)-(7) and, in
Sec. 5.2 the practical issue of setting the scale.

###  5.1 Effect of the knockout/shootout rules

Table 7 compares the ranking (of top-five teams) obtained using the FIFA
algorithm (first column) to the rating resulting from the modified algorithm
in which we a) eliminate the shootout rule (second column), b) eliminate the
knockout rule (third column), as well as c) elimination both rules (fourth
column). The differences, most notably the removal of Belgium from the first
place, are due to the different number of times the teams benefited from the
knockout rules (although the shootout rule for sure has an effect on the final
rating too).

Indeed, by analyzing the results of the games, we observed that in the
original ranking, Belgium (BEL) benefited four times from the knockout rule
for a total of 85 points (which would be lost without the rule (6)), Brazil
(BRA) and England (ENG) benefited twice for a total of 53 and 80 points,
respectively, while both France (FRA) and Italy (ITA) benefited only once,
gaining 14 points each.111111Of course, due to the temporal relationships,
eliminating the knockout/shootout rules is not the same as evaluating the
points (not lost in the original algorithm) and discarding them from the final
results.

What is important is that the points-preserving knockout rule ignores the
direct comparison between the teams. In fact, the games in which Belgium was
not penalized (for loosing in knockout stages) were played against France
(twice) and against Italy (twice as well). Thus, despite a direct evidence
indicating that France and Italy were able to beat Belgium, the knockout rule
preserved the points earned by Belgium in other games.

In fact, such a situation is not surprising and we indeed expect the teams
which compete for the top ranking spots to be also likely to make it to the
final stages of the important competitions (in case of the Belgium’s games:
World Cup 2018, Euro 2020, and UEFA Nations League 2021) and then play against
each other. While these games will provide direct comparison results, current
knockout rule will preserve the points of the losing team.

Whether this is fair and desirable may be debatable especially considering
that the knockout/shootout rules are not rooted in any formal modelling
principle, and most likely are introduced to compensate for the increased
value of $I_{c}$ in the advanced stages of competitions.

###  5.2 Scale adjustment

The scale is obviously irrelevant in the batch optimization and the on-line
update can also be written in the scale-invariant manner by dividing (18) by
$s$:

$\displaystyle\boldsymbol{\theta}^{\prime}_{t+1}$
$\displaystyle\leftarrow\boldsymbol{\theta}^{\prime}_{t}-K^{\prime}\xi_{c_{t}}\boldsymbol{x}_{t}g(z^{\prime}_{t};y_{t})$
(51) $\displaystyle z^{\prime}_{t}$ $\displaystyle=z_{t}/s$ (52)
$\displaystyle\boldsymbol{\theta}^{\prime}_{t}$
$\displaystyle=\boldsymbol{\theta}_{t}/s$ (53) $\displaystyle K^{\prime}$
$\displaystyle=K/s;$ (54)

in other words, for the same scale-invariant initialization
$\boldsymbol{\theta}^{\prime}_{0}$ and using the same step $K^{\prime}$ we
will obtain the same results $\boldsymbol{\theta}^{\prime}_{t}$.

However, in the FIFA ranking, a non-zero initialization
$\boldsymbol{\theta}_{0}$ was determined in advance and thus
$\boldsymbol{\theta}^{\prime}_{0}$ is not scale-invariant. Thus, given the
initialization at hand, the question is how to determine the scale? In
general, it is, of course, a difficult question but an insight may be gained
assuming that the initialization corresponds to the “optimal” solution, e.g.,
$\hat{\boldsymbol{\theta}}$ obtained in the batch optimization with a given
scale $s_{0}$.

It is easy to see that using $s>s_{0}$ will force the algorithm to
significantly change $\boldsymbol{\theta}_{t}$ (attainable with large values
of the adaptation step, $K$); the same will happen for $s<s_{0}$ because the
optimal estimates $\boldsymbol{\theta}_{t}$ will have to be scaled down.

Since scaling up/down of the skills changes their empirical moments we suggest
to choose the scale, $s$ in a moment-preserving manner. To this end we define
the empirical standard deviation of the skills

$\displaystyle\sigma_{t}=\sqrt{\|\hat{\boldsymbol{\theta}}_{t}-\overline{\hat{\boldsymbol{\theta}}}_{t}\|^{2}/M}$
(55)

where
$\overline{\hat{\boldsymbol{\theta}}}_{t}=(\sum_{m=1}^{M}\hat{\theta}_{t,m})/M$
is the empirical mean, and postulate that, at the initialization and at the
final step, we have $\sigma_{0}\approx\sigma_{T}$.

In fact, the initialization used by FIFA yields $\sigma_{0}=220$ and, after
running the FIFA algorithm we obtain $\sigma_{T}=250$, relatively close to the
initial value $\sigma_{0}$.

Changing the scale $s$, we will obtain different $\sigma_{T}$ so the idea is
to run the algorithms for different values of the scale $s$, e.g., as
multiples of $100$ and to choose the one which yields a standard deviation
$\sigma_{T}\approx\sigma_{0}$. In practice it has to be done using historical
data _before_ the new rating is deployed but in our case we could do it in the
hindsight.

In this manner we found $s=200$ to be suitable for the Davidson-Elo algorithm
(we obtained  $\sigma_{T}=219$ for the unweighted version and $\sigma_{T}=221$
for the MOV-weigted approach), and $s=300$ well suited for the Skellam’s
algorithms (where  $\sigma_{T}=225$ was obtained). This also indicates the the
scale $600$ was too large for the FIFA rating. This can be noted by comparing,
in Table 8 the result FIFA with $\xi_{c}=1$ to the results of the SG (with
$\eta=0$ and $\kappa=2$). Both are essentially the same algorithms (although
FIFA uses the shootout/knockout rules which have negligible impact on the
performance) and the only difference resides in the scale. Since the scale
$s=200$ in the Elo-Davidson algorithm corresponds to the scale $s=400$ in the
FIFA algorithm, the latter would perform better with the scale $s=400$. This
effect, however, appears only due to limited observation window we have at our
disposal and will vanish after a sufficiently large number of games.

algorithm | $\mathsf{LS}_{\textrm{opt}}$ | $K$ | $\eta$ | $\kappa$ | $\mathsf{ACC}~{}[\%]$
---|---|---|---|---|---
FIFA, $\xi_{c}$ from Table 1 | $0.951$ | $5$ | $0$ | $2$ | $50$
FIFA, $\xi_{c}=1$ | $0.933$ | $55$ | $0$ | $2$ | $52$
SG | $0.917$ | $35$ | $0$ | $2$ | $54$
$0.892$ | $35$ | $0.4$ | $2$ | $58$
$0.841$ | $35$ | $0.3$ | $0.9$ | $61$

a) Performance of the algorithms : FIFA ($s=600$) and Elo-Davidson model with
SG ($s=200$)  $\mathsf{LS}_{\textrm{opt}}$ $V$ $K$ $\eta$ $\kappa$ $\zeta_{0}$
$\zeta_{1}$ $\zeta_{2}$ $\zeta_{3}$ $\mathsf{ACC}~{}[\%]$ $0.841$ $1$ $35$
$0.3$ $0.9$ $1.0$ $0.9$ $\times$ $\times$ $61$ $0.838$ $2$ $35$ $0.3$ $0.9$
$1.0$ $0.6$ $1.3$ $\times$ $62$ $0.837$ $3$ $40$ $0.3$ $0.9$ $1.0$ $0.5$ $0.8$
$1.8$ $62$ b) SG with the MOV weighting, $s=200$ $\mathsf{LS}_{\textrm{opt}}$
$K$ $\eta$ $c$ $\mathsf{ACC}~{}[\%]$ $0.827$ $7.5$ $0.2$ $-0.1$ $62$ c) SG
implementing Skellam’s model for the MOV, $s=300$

Table 8: Parameters and performance of the on-line rating SG algorithms
obtained by minimizing the log-score (56) for a) Davidson model, b) MOV-
weighting strategy from Sec. 4.1, and c) the MOV-modelling strategy from Sec.
4.2.

###  5.3 Evaluation of the algorithms

To evaluate the SG algorithms for the models studied in the batch context, we
will use the log-score and the accuracy metrics defined for the half of the
games in the considered time period

$\displaystyle\mathsf{LS}$
$\displaystyle=\frac{2}{T}\sum_{t=T/2+1}^{T}\mathsf{m}^{\textrm{ls}}(z_{t},y_{t})$
(56) $\displaystyle\mathsf{ACC}$
$\displaystyle=\frac{2}{T}\sum_{t=T/2+1}^{T}\mathsf{m}^{\textrm{acc}}(z_{t},y_{t}).$
(57)

We consider the SG algorithm based on the Davidson model (Table 8a), the
Davidson model with the MOV-weighting (Table 8b), and Skellam’s model
algorithm (Table 8c).

In all cases, but in the original FIFA algorithm, we ignore the category-
weighting (i.e., we use $\xi_{c}=1$) because, as we have already shown, its
effect is negligible. This is clearly shown in the first part of Table 8a
where we see that using the FIFA weighting we obtain worse results than when
the weighting in ignored. This is essentially the same result as the one we
have shown in Table 2 but we repeat it here to show the log-score metric which
we could not calculate without first introducing the Davidson model underlying
the FIFA algorithm.

The results indicate that:

  * •

The most notable improvements are due to, in similar measures, two elements:
the introduction of the HFA coefficient, $\eta$ and the explicit use of the
Davidson model (and thus, the optimization of the coefficient $\kappa$).

  * •

Additional small, but still perceivable gains are obtained by introducing the
MOV-weighting, where from the lesson learnt in Sec. 4.1 we weight
independently the draws and the home/away wins.

  * •

The MOV-modelling using the Skellam’s distribution brings again a small
benefit.

We present in Table 9 the rating obtained for the top teams via new rating
algorithms. Of course, due to smaller scale we used, the skills have smaller
values and should not be compared directly to those from Table 7 but the
ranking is of interest, where the teams from the FIFA ranking are present
(FRA, BRA, BEL) but this time Argentine (ARG), which was on the sixth place in
the previous rankings, is now consistently on and above the top-third
position. We can also see that the differences between the rating values are
much less pronounced.

Davidson | MOV weights | Skellam’s model
---|---|---
FRA (1683.5) | FRA (1690.8) | BRA (1596.0)
BRA (1673.1) | BRA (1677.4) | ARG (1585.8)
ARG (1668.6) | ARG (1677.0) | BEL (1546.2)
BEL (1664.9) | BEL (1666.5) | POR (1541.2)
ITA (1657.7) | ITA (1665.6) | ESP (1540.7)
Table 9: Ranking of the top teams using the proposed algorithms.

##  6 Conclusions

In this work we analyzed the FIFA ranking using the methodology conventionally
used in the probabilistic modelling and inference. In the first step, we
identified the model relating the outcomes (games results) to the parameters
which have to be optimized (skills of the teams). More precisely, we have
shown that the FIFA algorithm can be formally derived as the stochastic
gradient (SG)optimization of the weighted maximum likelihood (ML)criterion in
the Davidson model (Davidson, 1970).

This first step allows us to define the performance metrics related to the
predictive performance of the algorithms we study. This is particularly
important in the case of the FIFA ranking algorithm because it does not model
the outcomes of the game but only explicitly specifies the expected score,
which is not sufficient to precisely evaluate the rating results. It also
allows us to apply the batch approach to rating and skills’ estimation. This
conventional machine learning strategy frees us from the considerations
related to the scale, initialization, or modeling of the skills’ dynamics.

Using the batch rating, we have shown that the game-category weighting is
negligible at best, and counterproductive at worst, which is the case of the
weighting used by the FIFA rating. This observation is interesting in its own
right because, while on one hand the concept of weighting is also used in the
rating literature, e.g., (Ley ., 2019), on the other, the literature does not
show any evidence that it is in any way beneficial and our findings
consistently indicate the contrary.

We next considered extensions of the algorithm by including the HFA and
optimizing the parameter responsible for the draws. These two elements seem to
be particularly important from the point of view of the performance of the
rating algorithm. While the HFA is a well-known element, already considered by
FIFA in FIFA (2007), the possibility of generalizing the Elo algorithm by
using the Davidson’s model, was only recently shown in Szczecinski Djebbi
(2020).

We also evaluated the possibility of using the margin of victory (MOV)given by
the goal differential, where we analyzed the weighting strategy and the
modelling based on the Skellam’s distribution. These two methods further
improve the results at the cost of higher complexity. Here, the formal
optimization strategy of the weighting parameters also yield interesting and
somewhat counter-intuitive results. Namely, we have shown that the games won
with small margin should have smaller weights than the tied games. This stands
in net contrast with the weighting strategies proposed before, e.g., in
Hvattum Arntzen (2010), Silver (2014), Kovalchik (2020) which use the
weighting with monotonically increasing functions of the margin.

Finally, we evaluated the heuristic shootout/knockout rules which are used in
the FIFA rating. Since their impact on the overall performance is small and
they may distort the relationship between the ratings of the strong teams
which often face each other in the final stages of the competitions, their
usefulness is questionable. In particular, eliminating the knockout rule would
strip Belgium from its first place position in the current FIFA ranking due to
multiple losses Belgium suffered against the current top teams (e.g., Italy,
France).

###  6.1 Recommendations

Given the analysis and the observations we made, if the FIFA rating was to be
changed, the following steps are recommended:

  1. 1.

Add the home-field advantage (HFA)parameter to the model because playing at
the home venue is a strong predictor for the victory. Not only this well-known
fact is already exploited in Women FIFA ranking but such a modification is
most likely the simplest and the least debatable element. In our view, it is
surprising that the current rating adopted in 2018 does not include the HFA.

  2. 2.

Use explicit model to relate the skills to the outcomes. Not only it would add
expressiveness providing the explicit predicted probability for each outcomes,
but it also improves the prediction results. Note that the rating algorithm
introduced recently by FIVB adopts such an approach and specifies the
probability for each of the game outcomes. In the context of the FIFA ranking,
the Davidson model we used in this work is an excellent candidate for that
purpose as it results in a natural generalization of the Elo algorithm,
preserving the legacy of the current algorithm.

  3. 3.

Remove the weighting of the games according to their assumed importance
because the data does not provide any evidence for their utility, or rather
provides the indication that the weighting in its current form is
counterproductive. If the concept of the game importance is of extra-
statistical nature (such as entertainment), it is preferable to diminish its
role, e.g., by shrinking the gap between the largest and the smaller values of
$\xi_{c}$ used.

  4. 4.

Remove the shootout and knockout rules which are not rooted in any sound
statistical principle.

As far as the knockout rule is concerned, while the intent to protect the
rating of the teams which manage to qualify to the knockout stage is clear, we
may argue that the penalty due to losing in the knockout game is aggravated by
the increased weighting of these games. Therefore, removing the weight, as we
postulate, would also eliminate the very reason to protect the teams’ points
with the knockout rule.

Regarding the shootout rule, a small frequency of events where it can be
applied and a marginal changes in the score imposed by the rule, make its
impact rather negligible. Its fairness is again debatable because there is
little evidence relating the skills of the teams to the outcome of the
shootout.

  5. 5.

If the rating was to consider the MOV, the simplest solution lies in weighting
the update step using the goal differential. On the other hand, the
modification based on the change of the model using the Skellam’s distribution
may cause numerical problems and the relatively small performance gains hardly
justify the added complexity.

On the other hand, the MOV may be added using alternative solutions similar to
those already considered in the Women’ teams FIFA ranking. Again, the latter
should be studied, e.g., using the methodology we used in this work and basing
the results on a formal probabilistic model.

##  Appendix A Approximate leave-one-out cross-validation

Our goal is to calculate in a simple manner the terms
$\boldsymbol{x}_{t}^{\top}\hat{\boldsymbol{\theta}}_{\backslash{t}}$ which
appear in the scoring function in (30) and in (31).

We start by approximating the maximum a posteriori (MAP) objective function
(33) using the Taylor series

$\displaystyle J_{\backslash{t}}(\boldsymbol{\theta})$
$\displaystyle=J(\boldsymbol{\theta})+\xi_{c_{t}}\log
L(\boldsymbol{x}_{t}^{\top}\boldsymbol{\theta}/s;y_{t})$ (A.1)
$\displaystyle\approx J(\hat{\boldsymbol{\theta}})+\xi_{c_{t}}\log
L(\boldsymbol{x}_{t}^{\top}\hat{\boldsymbol{\theta}}/s;y_{t})$
$\displaystyle\quad-\frac{\xi_{c_{t}}}{s}g_{t}\boldsymbol{x}_{t}^{\top}(\boldsymbol{\theta}-\hat{\boldsymbol{\theta}})+\frac{1}{2}(\boldsymbol{\theta}-\hat{\boldsymbol{\theta}})^{\top}\Big{[}\hat{\textnormal{{H}}}-\frac{\xi_{c_{t}}}{s^{2}}h_{t}\boldsymbol{x}_{t}\boldsymbol{x}_{t}^{\top}\Big{]}(\boldsymbol{\theta}-\hat{\boldsymbol{\theta}})$
(A.2)

where $g_{t}\equiv
g(\boldsymbol{x}^{\top}_{t}\hat{\boldsymbol{\theta}}/s;y_{t})$ is defined in
(23),

$\displaystyle\hat{\boldsymbol{\theta}}$
$\displaystyle=\mathop{\mathrm{argmin}}_{\boldsymbol{\theta}}J(\boldsymbol{\theta})$
(A.3)

is the optimal solution for all data, the Hessian at optimum is given by

$\displaystyle\hat{\textnormal{{H}}}$
$\displaystyle=\nabla^{2}_{\boldsymbol{\theta}}J(\boldsymbol{\theta})|_{\boldsymbol{\theta}=\hat{\boldsymbol{\theta}}}=\sum_{t\in\mathcal{T}}\frac{\xi_{c_{t}}}{s^{2}}h_{t}\boldsymbol{x}_{t}\boldsymbol{x}^{\top}_{t}+\frac{\alpha}{s^{2}}\boldsymbol{I},$
(A.4)

and we use second derivative $h_{t}\equiv
h(\boldsymbol{x}^{\top}_{t}\hat{\boldsymbol{\theta}}/s)$ with (Szczecinski
Tihon, 2021, Sec. IV)

$\displaystyle h(z)$
$\displaystyle=\frac{\,\mathrm{d}}{\,\mathrm{d}z}g(z;y)=\frac{\left(\ln
10\right)^{2}}{4}\frac{\kappa 10^{0.5(z+\eta b)}+4+\kappa 10^{-0.5(z+\eta
b)}}{(10^{0.5(z+\eta b)}+\kappa+10^{-0.5(z+\eta b)})^{2}}.$ (A.5)

By equating the gradient of (A.2) to zero, we find the approximate solution to
the optimization problem

$\displaystyle\hat{\boldsymbol{\theta}}_{\backslash{t}}$
$\displaystyle\approx\mathop{\mathrm{argmin}}_{\boldsymbol{\theta}}{J}_{\backslash{t}}(\boldsymbol{\theta})$
$\displaystyle=\hat{\boldsymbol{\theta}}+\frac{\xi_{c_{t}}g_{t}}{s}\Big{[}\hat{\textnormal{{H}}}-\frac{\xi_{c_{t}}}{s^{2}}h_{t}\boldsymbol{x}_{t}\boldsymbol{x}_{t}^{\top}\Big{]}^{-1}\boldsymbol{x}_{t}$
(A.6)

and the terms
$\boldsymbol{x}_{t}^{\top}\hat{\boldsymbol{\theta}}_{\backslash{t}},t\in\mathcal{T}$
which appear as arguments of the metrics (30) and (31) can be now calculated
efficiently for all $t\in\mathcal{T}$ once $\hat{\boldsymbol{\theta}}$ is
known (Rad Maleki, 2020)(Burn, 2020)

$\displaystyle\boldsymbol{x}_{t}^{\top}\hat{\boldsymbol{\theta}}_{\backslash{t}}$
$\displaystyle\approx\boldsymbol{x}_{t}^{\top}\hat{\boldsymbol{\theta}}+\frac{\xi_{c_{t}}g_{t}}{s}\boldsymbol{x}_{t}^{\top}\Big{[}\hat{\textnormal{{H}}}-\frac{\xi_{c_{t}}}{s^{2}}h_{t}\boldsymbol{x}_{t}\boldsymbol{x}_{t}^{\top}\Big{]}^{-1}\boldsymbol{x}_{t}$
(A.7)
$\displaystyle=\boldsymbol{x}_{t}^{\top}\hat{\boldsymbol{\theta}}+\frac{\xi_{c_{t}}g_{t}}{s}\boldsymbol{x}_{t}^{\top}\Big{[}\hat{\textnormal{{H}}}^{-1}+\frac{\xi_{c_{t}}h_{t}}{s^{2}-\xi_{c_{t}}h_{t}\boldsymbol{x}_{t}^{\top}\hat{\textnormal{{H}}}^{-1}\boldsymbol{x}_{t}}\hat{\textnormal{{H}}}^{-1}\boldsymbol{x}_{t}\boldsymbol{x}_{t}^{\top}\hat{\textnormal{{H}}}^{-1}\Big{]}\boldsymbol{x}_{t}$
(A.8)
$\displaystyle=\boldsymbol{x}_{t}^{\top}\hat{\boldsymbol{\theta}}+\frac{\xi_{c_{t}}g_{t}a_{t}s}{s^{2}-\xi_{c_{t}}h_{t}a_{t}},$
(A.9)

where
$a_{t}=\boldsymbol{x}_{t}^{\top}\hat{\textnormal{{H}}}^{-1}\boldsymbol{x}_{t}$
and to pass from (A.7) to (A.8) we used the matrix inversion lemma (Barber,
2012, Ch. A.1.8).

The advantage of this formulation is clear: instead of solving $T$ times the
optimization problem (32), we only need to solve once the optimization defined
in (A.3). In comparison with the latter, the remaining operations of the
inversion of the matrix ${\textnormal{{H}}}_{0}$ and the multiplication
required to calculate $a_{t},t\in\mathcal{T}$, have a very small complexity.

The identical approach may be used to apply the ALO to the problem (40) but,
we have to replace $\xi_{c_{t}}$ in (A.9) with $\xi_{c_{t}}\zeta_{v_{t}}$.

In order to apply the ALO to the problem (48) we need a second derivative of
(46) which is given by

$\displaystyle h(z)=\frac{\,\mathrm{d}}{\,\mathrm{d}z}g(z;d)$
$\displaystyle=\mathrm{e}^{c}(\mathrm{e}^{z+b\eta}+\mathrm{e}^{-z-b\eta}).$
(A.10)

## References

  * Amiguet (2010) Amiguet10_thesisAmiguet, M.  2010\.  Adaptively Weighted Maximum Likelihood Estimation of Discrete Distributions Adaptively weighted maximum likelihood estimation of discrete distributions .  Universtité de Lausanne.
  * Barber (2012) Barber12_BookBarber, D.  2012\.  Bayesian reasoning and Machine Learning Bayesian reasoning and machine learning.  Cambridge University Press.
  * Burn (2020) Burn20Burn, R.  2020\.  Optimizing Approximate Leave-one-out Cross-validation to Tune Hyperparameters Optimizing approximate leave-one-out cross-validation to tune hyperparameters.  ArXivabs/2011.10218.  http://arxiv.org/abs/2011.10218
  * Davidson (1970) Davidson70Davidson, RR.  1970\.  On Extending the Bradley-Terry Model to Accommodate Ties in Paired Comparison Experiments On extending the Bradley-Terry model to accommodate ties in paired comparison experiments.  Journal of the American Statistical Association65329317–328.  http://www.jstor.org/stable/2283595
  * Duda . (2001) Duda_bookDuda, RO., Hart, PE. Stork, DG.  2001\.  Pattern Classification Pattern classification.  John Wiley & Sons.
  * Egidi Torelli (2021) Egidi21Egidi, L. Torelli, N.  2021\.  Comparing Goal-Based and Result-Based Approaches in Modelling Football Outcomes Comparing goal-based and result-based approaches in modelling football outcomes.  Social Indicators Research1562801–813.  https://doi.org/10.1007/s11205-020-02293-z 10.1007/s11205-020-02293-z
  * Elo (1978) Elo78_BookElo, AE.  1978\.  The Rating of chessplayers, past and present The rating of chessplayers, past and present.  Arco Publishing Inc.
  * eloratings.net (2020) eloratings.neteloratings.net.  2020\.  World Football Elo Ratings. World football Elo ratings.  [Oct. 8, 2020]https://www.eloratings.net/
  * FIFA (2007) fifa_rating_WFIFA.  2007may.  FIFA/Coca-Cola Women’s World Ranking. FIFA/Coca-Cola women’s world ranking.  [November 12, 2021]https://digitalhub.fifa.com/m/3d9cb1decbbb2ac7/original/rxqyxdjhbs2qdtstluy6-pdf.pdf
  * FIFA (2018) fifa_ratingFIFA.  2018\.  Revision of the FIFA/Coca-Cola World Ranking. Revision of the FIFA/Coca-Cola world ranking.  [February 7, 2020]https://digitalhub.fifa.com/m/f99da4f73212220/original/edbm045h0udbwkqew35a-pdf.pdf
  * FIVB (2020) fivb_ratingFIVB.  2020\.  New Senior world rankings. New senior world rankings.  [December 6, 2021]https://www.fivb.com/en/volleyball/rankings
  * Football Rankings (2021) football_rankingsFootball Rankings.  2021\.  Football Rankings. Football rankings.  [October 28, 2021]http://www.football-rankings.info/
  * Gelman . (2014) Gelman14Gelman, A., Hwang, J. Vehtari, A.  20141101\.  Understanding predictive information criteria for Bayesian models Understanding predictive information criteria for Bayesian models.  Statistics and Computing246997–1016.  https://doi.org/10.1007/s11222-013-9416-2 10.1007/s11222-013-9416-2
  * Hastie . (2009) Hastie_bookHastie, T., Tibshirani, R. Friedman, J.  2009\.  The elements of Statistical Learning The elements of statistical learning.  Springer Series in Statiscs.
  * Hu Zidek (2001) Hu01Hu, F. Zidek, JV.  2001\.  The Relevance Weighted Likelihood With Applications The relevance weighted likelihood with applications.  SE. Ahmed N. Reid (), Empirical Bayes and Likelihood Inference Empirical bayes and likelihood inference ( 211–235).  New York, NYSpringer New York.  https://doi.org/10.1007/978-1-4613-0141-7_13 10.1007/978-1-4613-0141-7_13
  * Hvattum Arntzen (2010) Hvattum10Hvattum, LM. Arntzen, H.  2010\.  Using Elo ratings for match result prediction in association football Using Elo ratings for match result prediction in association football.  International Journal of Forecasting263460 - 470\.  http://www.sciencedirect.com/science/article/pii/S0169207009001708 Sports Forecasting  https://doi.org/10.1016/j.ijforecast.2009.10.002
  * Karlis Ntzoufras (2008) Karlis08Karlis, D. Ntzoufras, I.  200809\.  Bayesian modelling of football outcomes: using the Skellam’s distribution for the goal difference Bayesian modelling of football outcomes: using the Skellam’s distribution for the goal difference.  IMA Journal of Management Mathematics202133-145.  https://doi.org/10.1093/imaman/dpn026 10.1093/imaman/dpn026
  * Kovalchik (2020) Kovalchik20Kovalchik, S.  2020\.  Extension of the Elo rating system to margin of victory Extension of the Elo rating system to margin of victory.  International Journal of Forecasting.  http://www.sciencedirect.com/science/article/pii/S0169207020300157 https://doi.org/10.1016/j.ijforecast.2020.01.006
  * Lasek Gagolewski (2020) Lasek20Lasek, J. Gagolewski, M.  2020\.  Interpretable sports team rating models based on the gradient descent algorithm Interpretable sports team rating models based on the gradient descent algorithm.  International Journal of Forecasting.  http://www.sciencedirect.com/science/article/pii/S0169207020301849 https://doi.org/10.1016/j.ijforecast.2020.11.008
  * Lasek . (2013) Lasek13Lasek, J., Szlávik, Z. Bhulai, S.  2013\.  The predictive power of ranking systems in association football The predictive power of ranking systems in association football.  International Journal of Applied Pattern Recognition1127-46.  https://www.inderscienceonline.com/doi/abs/10.1504/IJAPR.2013.052339 PMID: 52339  10.1504/IJAPR.2013.052339
  * Ley . (2019) Ley19Ley, C., de Wiele, TV. Eetvelde, HV.  2019\.  Ranking soccer teams on the basis of their current strength: A comparison of maximum likelihood approaches Ranking soccer teams on the basis of their current strength: A comparison of maximum likelihood approaches.  Statistical Modelling19155-73.  https://doi.org/10.1177/1471082X18817650 10.1177/1471082X18817650
  * Maher (1982) Maher82Maher, MJ.  1982\.  Modelling association football scores Modelling association football scores.  Statistica Neerlandica363109-118.  https://onlinelibrary.wiley.com/doi/abs/10.1111/j.1467-9574.1982.tb00782.x 10.1111/j.1467-9574.1982.tb00782.x
  * Rad Maleki (2020) Rad20Rad, KR. Maleki, A.  2020\.  A scalable estimate of the out-of-sample prediction error via approximate leave-one-out cross-validation A scalable estimate of the out-of-sample prediction error via approximate leave-one-out cross-validation.  Journal of the Royal Statistical Society: Series B (Statistical Methodology)824965-996.  https://rss.onlinelibrary.wiley.com/doi/abs/10.1111/rssb.12374 https://doi.org/10.1111/rssb.12374
  * The Roon Ba (2021) roonbaThe Roon Ba. The Roon Ba.  2021\.  [Dec. 2, 2021]http://www.theroonba.com/
  * Silver (2014) Silver14Silver, N.  2014\.  Introducing NFL Elo ratings. Introducing NFL Elo ratings.  [July 1, 2020]https://fivethirtyeight.com/features/introducing-nfl-elo-ratings/
  * SoccerWay (2021) soccerwaySoccerWay. Soccerway.  2021\.  [Dec. 2, 2021]https://www.soccerway.com
  * Szczecinski (2020) Szczecinski20cSzczecinski, L.  2020\.  G-Elo: Generalization of the Elo algorithm by modelling the discretized Margin of Victory. G-Elo: Generalization of the Elo algorithm by modelling the discretized margin of victory.
  * Szczecinski Djebbi (2020) Szczecinski20Szczecinski, L. Djebbi, A.  2020\.  Understanding draws in Elo rating algorithm Understanding draws in Elo rating algorithm.  J. Quant. Anal. Sports163211-220.
  * Szczecinski Tihon (2021) Szczecinski21Szczecinski, L. Tihon, R.  2021\.  Simplified Kalman filter for online rating: one-fits-all approach. Simplified Kalman filter for online rating: one-fits-all approach.  http://arxiv.org/abs/2104.14012