# Seq2Seq Surrogates of Epidemic Models to Facilitate Bayesian Inference

Giovanni Charles  Imperial College London  Timothy M. Wolock  Imperial College
London  Peter Winskill  Imperial College London  Azra Ghani  Imperial College
London  Samir Bhatt  Seth Flaxman  Department of Computer Science, University
of Oxford

###### Abstract

Epidemic models are powerful tools in understanding infectious disease.
However, as they increase in size and complexity, they can quickly become
computationally intractable. Recent progress in modelling methodology has
shown that surrogate models can be used to emulate complex epidemic models
with a high-dimensional parameter space. We show that deep sequence-to-
sequence (seq2seq) models can serve as accurate surrogates for complex
epidemic models with sequence based model parameters, effectively replicating
seasonal and long-term transmission dynamics. Once trained, our surrogate can
predict scenarios a several thousand times faster than the original model,
making them ideal for policy exploration. We demonstrate that replacing a
traditional epidemic model with a learned simulator facilitates robust
Bayesian inference.

##  1 Introduction

Epidemic models are an essential tool in the study of infectious disease
dynamics [8, 26, 7]. By modelling new infections as a function of current
infections, they leverage the fundamental dynamics to predict epidemic
trajectories accurately with relatively few parameters. Recently, they have
risen to unprecedented public prominence through their use during the SARS-
CoV-2 (coronavirus) pandemic [3].

In most cases, we do not observe these dynamics or model parameters directly,
so we must infer them from what data we do measure (e.g. deaths, positive
tests, etc.). In order to produce a set of predictions, we need to infer the
unknown parameters in the model, a process that quickly becomes infeasible at
the scale needed for methods like Hamiltonian Monte Carlo (HMC) [15] for even
moderately complex models. Further, because the current state of the model is
a function of all past states, complex correlation structures can arise
between time-varying parameters, making sampling difficult.

[11] proposed an epidemiologically motivated Individual-Based Model (IBM) of
malaria transmission. It simulates infection at an individual-level based on
exposure to three species of infectious mosquito, immunity and a variety of
interventions. This model is invaluable for projecting the effects of policy
on malaria prevalence , for example, [25], projected the effects of a 44%
funding cut from a major global donor. Unfortunately, run times for a single
parameter set can take several minutes making it unsuitable for parameter
inference where large numbers of repeated function evaluations are required.
[10] later showed that a faster, reduced surrogate model can be used to infer
these parameters. This surrogate made simplifying assumptions, such as the
effects of seasonality and existing interventions being negligible. The
observed data had to be carefully curated to minimise those effects. In
contrast to a surrogate where a researcher arbitrarily chooses how to reduce
model complexity for inference, we will introduce a deep learning solution.

Recent work has demonstrated the utility of learned surrogate models reducing
the computational cost of simulation [4, 14, 17]. Data-driven methods could
potentially learn more expressive, surrogates which allow researchers to relax
their assumptions in parameter inference and use a wider range of data sets.
[4] proposed a functional regression method to map the functional space
between model parameters and outputs using a library of pre-computed outputs.
As the input space becomes large, more sophisticated methods have been
employed to explore the sub-spaces most useful for the parameter inference
task at hand, as reviewed by [14]. [17] showed that a Bayesian optimisation
workflow is effective in performing parameter inference on IBM outputs. They
proposed training a surrogate and simultaneously performing inference,
adaptively sampling parameter sets for training the surrogate which reduce the
surrogate’s predictive uncertainty in regions of the input space which are
proximal to the best-fit parameters.

However none of the above approaches have developed surrogates designed for
sequential inputs. Many epidemiological models use sequences as input, such as
transmissibility [16], intervention usage [10] and contact matrices [23], and
produce sequential outputs such as prevalence trajectories. This format makes
them amenable to a sequence-to-sequence (seq2seq) modelling approach. Advances
in seq2seq models have been very promising in machine translation [22] and
time-series forecasting [24], but it is yet to be seen if they can be applied
to the surrogate modelling of epidemiological models.

We hypothesised that a seq2seq based surrogate architecture could accurately
predict IBM outputs for a combination of scalar and sequential parameters.
Such a surrogate would make Bayesian inference feasible for long running
observational data sets.

###  1.1 Contributions

We created a surrogate for a well known malaria model, allowing researchers to
parameterise the surrogate for an arbitrary location and predict the malaria
prevalence within milliseconds. Previous approaches have either restricted the
input space to scalar parameters, or selected training samples for a specific
inference task. As far as we know, no other infectious disease model
surrogates are designed to take arbitrary sequence based parameters as input
at prediction time.

We propose a formulation for the surrogate in section 3. We outline a
surrogate model structure and training method in sections 3.3 and 3.4 where we
use synthetic parameter sets designed to explore the parameter space. We
evaluated the surrogate on two benchmark parameter sets, in section 4, to show
that it generalises to longer parameter sequences and parameter sets generated
from real world data.

We then performed parameter inference using a Hamiltonian Monte Carlo
algorithm on the learned surrogate. We estimated a parameter for a location in
Senegal, from an 18 year sequence of intervention statistics. This is the
first approach which would allow for inference using sequential parameters
over an arbitrary time horizon and/or location without re-training the
surrogate, which can reportedly take over a week per inference task, see [17].

##  2 Model of Malaria Transmission

We first introduce a malaria model to provide context for the formulation of
the surrogate model and to illustrate the complexity of an infectious disease
model.

We used an implementation of an Individual-Based Model (IBM) of malaria
originally published by [11]. The IBM includes human, mosquito and Insecticide
Treated Net (ITN) models. The human model simulates infectious bites,
immunity, and the probability of detection by diagnostic tests. This model is
used to produce malaria burden statistics and calculate the Force of Infection
towards Mosquitoes (FOIM), as explained in section 2.2. The mosquito model
simulates the size of the infectious and susceptible mosquito populations
given the FOIM and a time-series of rainfall. This model is used to calculate
an Entomological Inoculation Rate (EIR), explained in section 2.1, which
drives infections in the human model completing the cycle of infection between
humans and mosquitoes. The bed net model modifies this cycle by reducing the
number of bites in the protected human population and increasing mortality in
the mosquito population. These effects are explained in section 2.3.

###  2.1 Human Transmission

The human model, shown in figure 1, simulates five infection states:
susceptible (S), clinical disease (D), treated (T), asymptomatic infection
(A), and subpatent infection (U). Subpatent infections are those not detected
using routine malaria testing, instead requiring advanced diagnostic methods.
At each time step, Bernoulli trials are simulated to transition individuals
between these states at the parameterised rates.

Figure 1:  A human state transition diagram. Boxes represent human states and
arcs transitions between the states, simulated with Bernoulli trials.
$\Lambda$ represents the EIR, $f_{t}$ the probability of treatment, $b$ the
probability of infection and $\phi$ the probability of clinical disease.

The Entomological Inoculation Rate (denoted $\Lambda$), is calculated as the
number of infectious mosquito bites the population experiences over time.
However, each individual experiences $\Lambda$ differently. Older individuals
are bitten more often, due to their increased surface area. The age $a_{i}$ of
each individual $i$ is sampled from an exponential distribution with mean
$\mu_{a}$.

After an individual $i$ is bitten by an infectious mosquito, we simulate a
possible infection. $b_{i}(t)$ represents the probability of blood-stage
immunity failing, determining whether the individual will develop (at least)
an asymptomatic infection. $\phi_{i}(t)$ represents likewise probabilities for
clinical immunity and infection. For asymptomatic individuals,
$q_{i}(t,a_{i})$ represents the probability of detection by a routine
diagnostic test, taking into account the individual’s age. Each of these
immunity functions are complex estimations based on each individual’s previous
exposures to bites and infection.

The diagnostic prevalence of malaria in the model $\mbox{PfPr}(t,l,u)$ is the
proportion of individuals who would test positive in a routine diagnostic
test, as defined in equation 1. This statistic is bounded between lower and
upper age ranges $l$ and $u$, including all clinically diseased individuals
and asymptomatic individuals who would be detected.

$\begin{gathered}\mbox{PfPr}(t,l,u)=\frac{\sum_{i}\mbox{detected}(i,t)\times\mathds{1}_{l<a_{i}\leq
u}}{\sum_{i}\mathds{1}_{l<a_{i}\leq u}}\\\
\mbox{detected}(i,t)=\begin{cases}1&\mbox{state(i, t)}=D\\\
q_{i}(t,a)&\mbox{state}(i,t)=A\\\ 0&\mbox{otherwise}\end{cases}\end{gathered}$
(1)

###  2.2 Mosquito Modelling

Mosquitoes are modelled using a system of delay differential equations with
early larval (E), late larval (L), pupal (P), susceptible (Sm), incubating
(Em) and infectious (Im) compartments. Early (E) and Late (L) larval
compartments represent the number of larvae soon after hatching. The size of
these two compartments are dependent on a carrying capacity $K$ which can vary
with rainfall $R(t)$, see equation 3. The surviving larvae develop into Pupae
(P) and half of the surviving pupae will develop into susceptible female
mosquitoes.

$\displaystyle K^{v}=K^{v}_{0}\frac{R(t)}{\bar{R}}$ (2) $\displaystyle
R(t)=g_{0}+\sum_{i=1}^{3}g_{i}cos(2\pi ti)+h_{i}sin(2\pi ti)$ (3)

Susceptible female mosquitoes (Sm) are infected based on FOIM. FOIM is a
driver of infection in the mosquito population, calculated by the number of
bites each mosquito would have made on an infectious human. Infected
mosquitoes will incubate in the (Em) state for a fixed period before
themselves becoming infectious. If they survive, they will transition to the
Infected (Im) state where they remain until death. The number of new E state
larvae introduced each time step is proportional to the total number of adult
female mosquitoes in the simulation.

The mosquito model is replicated for each species of mosquito in the
simulation $v$. Users can parameterise the proportion of adult female
mosquitoes in each species $\kappa^{v}$ by setting proportional values of the
carrying capacity parameter, $K^{v}_{0}$. It is often convenient to
parameterise $K^{v}_{0}$ using a baseline EIR parameter, $\Lambda_{0}$, or the
expected EIR before interventions. [9] explains in detail how to calculate the
carrying capacity from baseline EIR.

###  2.3 ITN Modelling

Insecticide treated bednets (ITNs) are an effective malaria intervention which
protect individuals while they are in bed. They are treated with an
insecticide which repels or kills a mosquito when they attempt to feed on a
human while in bed. Mosquitoes which are repelled will take longer to feed,
indirectly affecting their mortality and reproduction rates. These responses,
for various species and insecticides, have been examined in experimental hut
trials, [18]. On the human side, EIR is reduced, since mosquitoes which have
been repelled or killed before feeding no longer contribute towards infectious
bites. We assume that bed nets are retained for an average of 5 years before
being thrown away. The model can be parameterised with a sequence, $\nu(t)$,
describing the proportion of the population who receive an ITN each year.

###  2.4 Model Fitting Considerations

Typically, malaria researchers would like to infer the unknown $\Lambda_{0}$
from observed prevalence data $\mbox{PfPr}(t,l,u)$, given known parameters
about the transmission setting, such as $\mu_{a},\kappa^{v},\nu(t),g,h...$.
However, simulating each individual, the mosquito population and ITN effects
for several years can take several minutes, making parameter inference
challenging. We also note that the individuals simulated by this model are
synthetic and not measured directly by any real-world data source. Any fitting
or inference process for this model will aggregate over individuals to predict
population-level averages which can be compared to data.

In the next section we will formulate a model to predict these population-
level aggregates directly from the IBM parameters, in order to eliminate a
large amount of computation that is, for the purposes of fitting the model to
data, unnecessary.

##  3 Surrogate Model

Let $f:\mathcal{X}\to\mathcal{Y}$ be a model to be emulated. Assume that
$\mathcal{X}\subset\mathbb{R}^{k}$ can be bounded to some set of plausible
values. Our goal is to train a model
$\tilde{f}_{\theta}:\mathcal{X}\to\mathcal{Y}$ such that
$|\tilde{f}_{\theta}(\mathbf{x})-f(\mathbf{x})|<\epsilon$ for some suitably
small $\epsilon>0$. We refer to $\tilde{f}_{\theta}$ as a surrogate model. Let
$\mathcal{L}(\mathbf{x},\mathbf{y},\tilde{f}_{\theta})$ be the loss of the
surrogate model with respect to $\mathbf{y}$ (where
$\mathbf{y}=f(\mathbf{x})$).

The goal of surrogate modelling is to find a parameter set that replicates the
original model accurately over the entirety of $\mathcal{X}$. That is, we want
to find $\hat{\theta}$:

$\hat{\theta}=\underset{\theta}{\mathrm{arg\;min}}\int_{\mathcal{X}}\mathcal{L}(\mathbf{x},\mathbf{y},\tilde{f}_{\theta})p(\mathbf{x})d\mathbf{x}.$
(4)

We have included a prior distribution over $\mathcal{X}$, but it is not
strictly necessarily. Here, $p(\mathbf{x})$ serves as a measure of the
relative importance of emulating $f(\mathbf{x})$ accurately. If $p$ is
uniform, then we are ascribing equal importance to every possible parameter
set.

###  3.1 Input Domain

We can translate the malaria IBM framing described in section 2 by first
defining $\mathcal{X}$.

To parameterise the surrogate for a particular location, we set $\kappa^{v}$,
$\mu_{a}$, $g_{i},h_{i}$ and $\nu(t)$ parameters to capture mosquito
composition, demography, seasonality and ITN usage respectively. We defined a
range for $\mu_{a}$ as 14.8 to 55.4 years based on the maximum and minimum
average age statistics per country from World Population Review [28]. Ranges
for Fourier coefficients $g_{i},h_{i}$ are set between -10 and 10, to include
the largest values found by fitting to rainfall data in sub-national
administrative units in sub-Saharan Africa. Since $\kappa^{v}$ are relative
proportions of three mosquito species, they are only valid in the 3-simplexes
in $\mathbb{R}^{3}$.

A realistic range of $\nu(t)$ is between 0 and 0.8 for each year, as it
becomes operationally difficult to distribute to ITNs to all individuals. We
chose not to constrain the relationship between $\nu(t)$ and $\nu(t+\delta)$
as it could limit the the effectiveness of the surrogate in exploratory
settings. However, parameterisations of $\nu(t)$ based on real world
statistics are likely to have a gradual increase from 0 and plateau. In
summary, $\nu(t)$ is defined as a bounded step function with $n\geq 0$
annually spaced knots.

One of the main unknowns for a each location is, $\Lambda_{0}$. An accurate
estimate of this parameter would allow us to simulate a population of humans
and mosquitoes which would recreate malaria transmission dynamics before any
interventions were implemented. Previous estimates for countries in sub-
Saharan Africa have ranged between 0.05 and 418 [10], and model outputs
suggest that higher values of $\Lambda_{0}$ have diminishing marginal effects
on PfPr. Therefore we set a range for $\Lambda_{0}$ between 0 (exclusive) and
500.

To express $\mathcal{X}$ as a sequence, we formulated a sequence of $n$
vectors for each year, $t$, with time-invariant parameters
($\kappa^{v},\mu_{a},g_{i},h_{i},\Lambda_{0})$ repeated at each time point,
and the $\nu(t)$ selected for the corresponding time point. Each vector
belongs to a subset of $\mathbb{R}^{13}$ bounded as described above.

###  3.2 Output Range

Let $(1,\ldots,T)$ be the sequence of years at which we wish to predict the
prevalence of malaria. For any $\mathbf{x}\in\mathcal{X}$, let $f$ represent
running the IBM with parameters $\mathbf{x}$ for $T$ years and evaluating
$\mbox{PfPr}(d,l,u)$ for each day $d$ in each year $t$:

$f(\mathbf{x})=\left\\{\left\\{\mbox{PfPr}(d,u,l)\right\\}_{d=1}^{365}\right\\}_{t=1}^{T}$
(5)

Note that we have augmented the notation for the IBM run to include parameter
sets. With this notation, $\mathcal{Y}$ is easy to define as the image of $f$:

$\mathcal{Y}=\\{f(\mathbf{x}):\mathbf{x}\in\mathcal{X}\\},$ (6)

where $\mathcal{X}$ is defined as in the previous section. From inspecting
equation 1, we get $\mathcal{Y}$ as a space of $T$-length sequences of vectors
$\mathbb{R}^{365}$ bounded between 0 and 1.

###  3.3 Surrogate Model Structure

With these definitions of $\mathcal{X}$ and $\mathcal{Y}$, the surrogate
modelling task can be reduced to a sequence-to-sequence prediction problem.
Sequence-to-sequence prediction is supported by a number of models, however we
chose a recurrent neural network (RNN) for simplicity. The RNN should be able
to exploit the relationships between time points to generalise for unseen
parameter sets and longer time horizons. We used Mean Squared Error (MSE) as
the loss function $\mathcal{L}$. We implemented a grid search select hyper-
parameters which minimised MSE on the validation set, table 1 lists the
considered parameters.

Hyper-parameter | Values
---|---
Optimiser |  Adam, RMSProp
Optimiser Loss | MSE, Log cosh
Batch size | 1000, 100, 50
Dropout |  0, .1
RNN Layer |  LSTM, GRU
Table 1: Hyper-parameters used in the grid search. The selected parameter is
emboldened

Our RNN took each vector from the $\mathcal{X}$ sequence as input to its
recurrent cells and output the corresponding vector from $\mathcal{Y}$. The
parameters for $g_{i}$ and $h_{i}$ were expanded to full rainfall data $R(t)$
using equation 3. The resulting inputs were standardised, removing the mean
and scaling to unit variance. The final input layer had 383 recurrent cells,
feeding 365 latent recurrent cells followed by a dense output layer of 365
cells.

###  3.4 Sampling Training Data

Training a surrogate for the malaria IBM requires a large library of
simulations covering a wide area of input space. We created training and
validation data, with an 80/20 split for the surrogate model using a sequence
length of $n=16$. This would allow for the surrogate to learn the long term
effects of ITNs on prevalence. We sampled 10,000 values from a 28 dimensional
latin-hypercube, 16 for $\nu(t)$ and 12 for the remaining time-invariant
parameters, and transformed the margins of each dimension to fit the bounds
listed in section 3.1. The IBM was run at the sampled points, setting aside a
warm-up period before activating the ITN parameters to allow the simulation to
reach an equilibrium. We set a large population size for the IBM runs to
minimise the stochasticity in PfPr trajectories and treat them as
deterministic outputs.

The simulation output was aggregated to calculate the prevalence of malaria
among children aged 5 to 59 months for each year after a warm-up period. The
age range of 5 to 59 months old was selected to align with statistics commonly
found in observational data.

##  4 Surrogate Evaluation

We evaluated the surrogate on two benchmark parameter sets. long contains
5,000 samples generated similarly to the training set except for sequences
with $n=30$, testing the surrogate’s ability to generalise to longer time
horizons. historic contains 5,000 samples where $\nu(t)$ was sampled from
historic statistics for a randomly selected country in sub-Saharan Africa,
published by [2], testing the performance on parameters representative of
global malaria modelling studies. After 100 epochs of the Adam optimiser, the
MSE across the 2,000 out-of-sample simulations was $6.8\times 10^{-3}$ on un-
normalised prevalence values. The prediction error for these benchmarks is
shown in figure 2.

Figure 2: True versus predicted prevalence of malaria in children aged 5-59
months (PfPr5_59) for the long (left) and historic (right) datasets where the
true estimates come from the IBM model and the predicted estimates come from
the surrogate. Each point represents the prevalence for one time step (a year)
of one sample in the corresponding benchmark. Mean Squared Error (mse) was
very low, as shown in the titles of the plots.

We compared execution times for the original IBM, surrogate training, single
prediction and batch prediction. The IBM executions were run on a cluster in
batches of 10 on 1,000 dedicated cores with a 2GB memory reservation. The
surrogate was trained on an Nvidia A1000 GPU, which took 229 $\pm$ 34 seconds
(n = 3). Both single and batch predictions were made on an 8-core Intel CPU.
While the IBM required just under 24 hours to run 10,000 samples, the batch
surrogate required just over 3 seconds (approximately a 25,000-fold
improvement). The timing statistics for an each run are listed in table 2.

Task | Seconds per run | Std
---|---|---
IBM run | $1.76\times 10^{3}$ | $7.69\times 10^{2}$
Prediction (single) | $1.74\times 10^{-3}$ | $2.84\times 10^{-8}$
Prediction (batch) | $3.74\times 10^{-4}$ | $6.12\times 10^{-2}$
Table 2: The timings per simulation for 10,000 simulations on the training and
validation parameter sets

##  5 Parameter Inference

We used Tensorflow Probability [20] to fit our model to observational data
from the Demographic and Health Survey program [6]. The DHS responses included
geo-location, age, and outcome of Rapid Diagnostic or Microscopy tests for
respondents in Kolda, Senegal between 2008 an 2018. We included responses for
respondents between 5 and 59 months old and weighted responses by population
density data from World Bank data [27] to get a general picture of parasite
prevalence over time in the region, $D_{\mbox{PfPr}5-59}$.

To estimate a value of $\Lambda_{0}$ which best explains
$D_{\mbox{PfPr}5-59}$, we used a binomial likelihood function, see equation 7.
We compared each participant’s infection status $J^{i}$, an indicator
variable, to the average predicted prevalence for the month in which they were
observed, $\mbox{PfPr}(t^{i})$. Our surrogate, parameterised with
$\nu(t),g_{i},h_{i},\kappa^{v}$ values for Kolda, and aggregated at month
$t^{i}$ provides quick estimates for $\mbox{PfPr}(t^{i})$. We set a uniform
prior for $\Lambda_{0}$ between 0 and 500.

$L(\lambda^{*}|D_{\mbox{PfPr}5-59})=\prod\mbox{PfPr}(t^{i})^{J^{i}}(1-\mbox{PfPr}(t^{i}))^{1-J^{i}}$
(7)

We performed Markov Chain Monte Carlo using the HMC algorithm [15] with 10
parallel chains, 2,000 steps, 1,000 of which were discarded as warm-up, to
sample from the resulting posterior. This inference took just over 11 minutes
to execute on an 8-core Intel CPU.

The original IBM was then run with the mean estimate for comparison, see the
posterior predictive check in figure 3 and resulting fit in figure 4. For this
fit, there is a noticeable error between the surrogate and the original model.
The seasonal cycle in prevalence is delayed and flattened, while the average
prevalence is underestimated. However, the relative effects of ITNs appear
similar between the models.

Figure 3: A posterior predictive check plotting $\Lambda_{0}$ samples from 10
MCMC chains. A previous rough estimate of baseline EIR in three Senegalese
populations was 1, 20 and 200 (infectious bites per person per year) for the
cities of Dakar, Ndiop and Dielmo [21]. This estimate of $23.2\pm 13.7$ would
place Kolda’s baseline transmission alongside Ndiop. Figure 4: A plot of the
final fit including the surrogate trajectory for the mean EIR (in orange). The
IBM was then calibrated to the mean EIR using a root finding method and run to
simulate the ITN distributions for Kolda, the prevalence for that run is
plotted in green. Observations have been plotted in blue with their 95%
confidence interval, assuming a Binomial distribution. The observational data
included samples from low intervention periods (2008) as well as high
intervention periods (2013-2017) and was aggregated on a monthly time scale.

##  6 Discussion

We have demonstrated that it is feasible to perform inference on a highly
detailed epidemic model with complex temporal dynamics. Seq2seq based
surrogates can accurately predict trajectories for models with large parameter
spaces and generalise transmission dynamics to very long time horizons. This
opens up a wider range of data, including longitudinal, from less controlled
studies for use in parameter inference in an epidemiological setting.

However, the inference in the previous section should not be used to inform
real policy. Simplifying assumptions have been made which ignore the effects
of non-ITN interventions and different diagnostic methods. Also, the accuracy
of the surrogate varies for different types of model runs. Further work needs
to be done to add detail to the surrogate model and ensure a high accuracy for
the parameter space which is used for inference.

###  6.1 A Surrogate for General Analysis

The surrogate could be used as an exploratory tool, allowing researchers to
begin modelling malaria transmission for arbitrary locations in seconds.
Selective approaches, such as [14, 17], had designed their training set for a
particular inference task, thus more training would be required to guarantee
their published performance for each use case. In the case of [17], training
can take between 7 and 12 days.

The benchmarks showed that the surrogate could emulate the IBM for the vast
majority of simulation scenarios. 99% and 97% of predictions had an MSE of
less than $1\times 10^{-2}$ for the long and historic benchmarks respectively.
On inspection, we could see that short-term seasonal and long-term ITN
transmission dynamics could be accurately re-created (see Figure 5).

Figure 5: 16 trajectories with the lowest MSE from the historic benchmark. The
blue and orange trajectories show daily prevalence in simulation over 16 years
for the IBM and surrogate models respectively.

Despite the benchmarks showing a small risk in surrogate prediction error,
many predictions significantly deviated from the IBM outputs. Manual
inspection reveals that the worst predictions for the historic benchmark were
in scenarios with low or zero prevalence (see Figure 6). The performance of
the surrogate would likely be improved by extending the training set with low
transmission parameter sets, i.e. samples with smaller values of $\Lambda$.

A more general approach to managing the risk in surrogate prediction error
could be uncertainty estimation. Unlike the approach in [17], RNNs provide
point estimates of the IBM outputs. Recent advances have shown that RNNs can
be adapted to output probability densities, [13, 12]. Such outputs could be
used at prediction time to inform the researcher of the trust-worthiness of a
prediction. They could also be used at training time to find samples which
would reduce surrogate predictive uncertainty, [1].

Figure 6: 16 cases with the highest MSE from the historic benchmark. The blue
and orange trajectories show daily prevalence in simulation over 16 years for
the IBM and surrogate models respectively.

Finally, sharing a surrogate model is more energy efficient and equitable than
an IBM. Low and middle income countries are disproportionately affected by
infectious disease and are less able to finance high performance computing.
Surrogate predictions require a small fraction of the IBM’s computational
requirements making analyses more feasible in the settings where they are most
needed.

###  6.2 Generality of the Workflow

Despite the overall input domain being large, the vectors for each time point
are low-dimensional. Over 50 malaria parameters were fixed to values sources
from expert opinion and previous parameter inference studies [19, 18, 10]. If
incorporated, it is unlikely that a latin-hypercube sampling strategy would
produce the samples sufficient to learn the function $f$ over the larger input
space. Developing a strategies to more intelligently sample the input space
without sacrificing generality poses a challenging research problem.

The sequence structure is also limited in granularity. This surrogate has
fixed time-step intervals to a year. However many infectious diseases,
including malaria, are seasonal and so there is value in exploring the timing
of interventions on a monthly or daily timescale. Reducing the time-step
interval would increase the size of training sequences. For example, training
this malaria surrogate with the daily time-step would increase the training
sequences from $n=16$ to $n=5840$. More advanced seq2seq architectures could
be considered here, such as those in [5, 22], taking care to avoid the
vanishing/exploding gradient problems.

##  7 Conclusion

We believe that the model presented here can be a valuable tool for guiding
future malaria modelling, and serves as a proof of concept for the utility of
surrogate modelling in infectious disease epidemiology. Given the limitations
we have explored, future work will focus on improving the surrogate’s fit and
generalising the surrogate structure to adapt it to more varied prediction
tasks.

## Acknowledgements

This research was funded in whole, or in part, by the Wellcome Trust [Grant
number 220900/Z/20/Z]. For the purpose of open access, the author has applied
a CC BY public copyright licence to any Author Accepted Manuscript version
arising from this submission.

S.B. acknowledges support from the MRC Centre for Global Infectious Disease
Analysis (MR/R015600/1), jointly funded by the UK Medical Research Council
(MRC) and the UK Foreign, Commonwealth & Development Office (FCDO), under the
MRC/FCDO Concordat agreement, and also part of the EDCTP2 programme supported
by the European Union. SB is funded by the National Institute for Health
Research (NIHR) Health Protection Research Unit in Modelling and Health
Economics, a partnership between the UK Health Security Agency, Imperial
College London and LSHTM (grant code NIHR200908). Disclaimer: “The views
expressed are those of the author(s) and not necessarily those of the NIHR, UK
Health Security Agency or the Department of Health and Social Care.” S.B.
acknowledges support from the Novo Nordisk Foundation via The Novo Nordisk
Young Investigator Award (NNF20OC0059309). SB acknowledges support from the
Danish National Research Foundation via a chair grant. S.B. acknowledges
support from The Eric and Wendy Schmidt Fund For Strategic Innovation via the
Schmidt Polymath Award (G-22-63345) which also supports GDC. SF acknowledges
the EPSRC (EP/V002910/2).

## References

  * [1] Tao Bai et al.  “Recent Advances in Adversarial Training for Adversarial Robustness” arXiv:2102.01356 [cs]  arXiv, 2021  URL: http://arxiv.org/abs/2102.01356
  * [2] S. Bhatt et al.  “The effect of malaria control on Plasmodium falciparum in Africa between 2000 and 2015”  In _Nature_ 526.7572, 2015, pp. 207–211  DOI: 10.1038/nature15535
  * [3] Paul Birrell et al.  “Real-time nowcasting and forecasting of COVID-19 dynamics in England: the first wave” Publisher: Royal Society  In _Philosophical Transactions of the Royal Society B: Biological Sciences_ 376.1829, 2021, pp. 20200279  DOI: 10.1098/rstb.2020.0279
  * [4] Ewan Cameron et al.  “Defining the relationship between infection prevalence and clinical incidence of Plasmodium falciparum malaria” Number: 1 Publisher: Nature Publishing Group  In _Nature Communications_ 6.1, 2015, pp. 8170  DOI: 10.1038/ncomms9170
  * [5] Kyunghyun Cho, Bart Merrienboer, Dzmitry Bahdanau and Yoshua Bengio  “On the Properties of Neural Machine Translation: Encoder-Decoder Approaches” arXiv: 1409.1259  In _arXiv:1409.1259 [cs, stat]_ , 2014  URL: http://arxiv.org/abs/1409.1259
  * [6] DHS  “The DHS Program - Quality information to plan, monitor and improve population, health, and nutrition programs” Accessed: 2022-05-12, https://dhsprogram.com/, 2022  URL: https://dhsprogram.com/
  * [7] Bimandra A Djaafara et al.  “A Quantitative Framework for Defining the End of an Infectious Disease Outbreak: Application to Ebola Virus Disease”  In _American Journal of Epidemiology_ 190.4, 2021, pp. 642–651  DOI: 10.1093/aje/kwaa212
  * [8] Jeffrey W. Eaton et al.  “The Estimation and Projection Package Age-Sex Model and the r-hybrid model: new tools for estimating HIV incidence trends in sub-Saharan Africa”  In _AIDS (London, England)_ 33.Suppl 3, 2019, pp. S235–S244  DOI: 10.1097/QAD.0000000000002437
  * [9] Jamie T. Griffin  “Is a reproduction number of one a threshold for Plasmodium falciparum malaria elimination?”  In _Malaria Journal_ 15.1, 2016, pp. 389  DOI: 10.1186/s12936-016-1437-9
  * [10] Jamie T. Griffin, Neil M. Ferguson and Azra C. Ghani  “Estimates of the changing age-burden of Plasmodium falciparum malaria disease in sub-Saharan Africa” Number: 1 Publisher: Nature Publishing Group  In _Nature Communications_ 5.1, 2014, pp. 3136  DOI: 10.1038/ncomms4136
  * [11] Jamie T. Griffin et al.  “Reducing Plasmodium falciparum Malaria Transmission in Africa: A Model-Based Evaluation of Intervention Strategies” Publisher: Public Library of Science  In _PLOS Medicine_ 7.8, 2010, pp. e1000324  DOI: 10.1371/journal.pmed.1000324
  * [12] Balaji Lakshminarayanan, Alexander Pritzel and Charles Blundell  “Simple and Scalable Predictive Uncertainty Estimation using Deep Ensembles” arXiv:1612.01474 [cs, stat]  arXiv, 2017  URL: http://arxiv.org/abs/1612.01474
  * [13] David J.. MacKay  “A practical Bayesian framework for backpropagation networks”  In _Neural Computation_ 4.3, 1992, pp. 448–472  DOI: 10.1162/neco.1992.4.3.448
  * [14] Trevelyan J. McKinley et al.  “Approximate Bayesian Computation and Simulation-Based Inference for Complex Stochastic Epidemic Models” Publisher: Institute of Mathematical Statistics  In _Statistical Science_ 33.1, 2018, pp. 4–18  URL: https://www.jstor.org/stable/26770975
  * [15] Radford M. Neal  “MCMC using Hamiltonian dynamics” arXiv: 1206.1901  In _arXiv:1206.1901 [physics, stat]_ , 2011  DOI: 10.1201/b10905
  * [16] Iman Rahimi, Fang Chen and Amir H. Gandomi  “A review on COVID-19 forecasting models”  In _Neural Computing and Applications_ , 2021  DOI: 10.1007/s00521-020-05626-8
  * [17] Theresa Reiker et al.  “Emulator-based Bayesian optimization for efficient multi-objective calibration of an individual-based model of malaria” Number: 1 Publisher: Nature Publishing Group  In _Nature Communications_ 12.1, 2021, pp. 7212  DOI: 10.1038/s41467-021-27486-z
  * [18] Ellie Sherrard-Smith et al.  “Optimising the deployment of vector control tools against malaria: a data-informed modelling study” Publisher: Elsevier  In _The Lancet Planetary Health_ 6.2, 2022, pp. e100–e109  DOI: 10.1016/S2542-5196(21)00296-5
  * [19] Marianne E. Sinka et al.  “A global map of dominant malaria vectors”  In _Parasites & Vectors_ 5.1, 2012, pp. 69  DOI: 10.1186/1756-3305-5-69
  * [20] TensorflowProbability  “TensorFlow Probability” Accessed: 2022-07-29  In _TensorFlow_ , https://www.tensorflow.org/probability/overview, 2022  URL: https://www.tensorflow.org/probability/overview
  * [21] J-F. Trape and C. Rogier  “Combating malaria morbidity and mortality by reducing transmission”  In _Parasitology Today_ 12.6, 1996, pp. 236–240  DOI: 10.1016/0169-4758(96)10015-6
  * [22] Ashish Vaswani et al.  “Attention Is All You Need” arXiv: 1706.03762  In _arXiv:1706.03762 [cs]_ , 2017  URL: http://arxiv.org/abs/1706.03762
  * [23] Oliver J Watson et al.  “Global impact of the first year of COVID-19 vaccination: a mathematical modelling study”  In _The Lancet Infectious Diseases_ , 2022  DOI: 10.1016/S1473-3099(22)00320-6
  * [24] Qingsong Wen et al.  “Transformers in Time Series: A Survey” arXiv: 2202.07125  In _arXiv:2202.07125 [cs, eess, stat]_ , 2022  URL: http://arxiv.org/abs/2202.07125
  * [25] Peter Winskill et al.  “The US President’s Malaria Initiative, Plasmodium falciparum transmission and mortality: A modelling study” Publisher: Public Library of Science  In _PLOS Medicine_ 14.11, 2017, pp. e1002448  DOI: 10.1371/journal.pmed.1002448
  * [26] World Health Organization  “Global Technical Strategy For Malaria - 2016 - 2030” Section: liv, 263 p.  Geneva: World Health Organization, 2021  URL: https://apps.who.int/iris/handle/10665/350147
  * [27] WorldBank  “World Bank Open Data — Data” Accessed: 2022-07-29, https://data.worldbank.org/, 2022  URL: https://data.worldbank.org/
  * [28] WorldPopulationReview  “2022 World Population by Country” Accessed: 2022-07-29, https://worldpopulationreview.com/, 2022  URL: https://worldpopulationreview.com/


11institutetext: pmOne Group
Paderborn, Germany 22institutetext: Institute of Informatics, University of
Munich (LMU)
33institutetext: Munich Center for Machine Learning (MCML)
Munich, Germany

Authors’ Instructions

# Weighting by Tying: A New Approach to Weighted Rank Correlation

Sascha Henzgen and Eyke Hüllermeier 112233

###### Abstract

Measures of rank correlation are commonly used in statistics to capture the
degree of concordance between two orderings of the same set of items. Standard
measures like Kendall’s tau and Spearman’s rho coefficient put equal emphasis
on each position of a ranking. Yet, motivated by applications in which some of
the positions (typically those on the top) are more important than others, a
few weighted variants of these measures have been proposed. Most of these
generalizations fail to meet desirable formal properties, however. Besides,
they are often quite inflexible in the sense of committing to a fixed weighing
scheme. In this paper, we propose a weighted rank correlation measure on the
basis of fuzzy order relations. Our measure, called _scaled gamma_ , is
related to Goodman and Kruskal’s gamma rank correlation. It is parametrized by
a fuzzy equivalence relation on the rank positions, which in turn is specified
conveniently by a so-called scaling function. This approach combines soundness
with flexibility: it has a sound formal foundation and allows for weighing
rank positions in a flexible way.

##  1 Introduction

Rank correlation measures such as Kendall’s tau [11] and Spearman’s rho [19],
which have originally been developed in non-parametric statistics, are used
extensively in various fields of application, ranging from bioinformatics [1]
to information retrieval [20]. In contrast to numerical correlation measures
such as Pearson correlation, rank correlation measures are only based on the
ordering of the observed values of a variable. Thus, measures of this kind are
not limited to numerical variables but can also be applied to non-numerical
variables with an ordered domain (i.e., measured on an ordinal scale) and, of
course, to rankings (permutations) directly.

In many applications, such as Internet search engines, one is not equally
interested in all parts of a ranking. Instead, the top positions of a ranking
(e.g., the first 10 or 50 web sites listed) are typically considered more
important than the middle part and the bottom. Standard rank correlation
measures, however, put equal emphasis on all positions. Therefore, they cannot
distinguish disagreements in different parts of a ranking. This is why
_weighted_ variants have been proposed for some correlation measures, as well
as alternative measures specifically focusing on the top of a ranking [6, 8,
10, 15, 20]. Most of these generalizations fail to meet desirable formal
properties, however. Besides, they are often quite inflexible in the sense of
committing to a fixed weighing scheme.

In this paper, we develop a general framework for designing weighted rank
correlation measures based on the notion of _fuzzy order relation_ , and use
this framework to generalize Goodman and Kruskal’s gamma coefficient [9]. Our
approach has a sound formal foundation and allows for weighing rank positions
in a flexible way. In particular, it is not limited to monotone weighing
schemes that emphasize the top in comparison to the rest of a ranking. The key
ingredients of our approach, to be detailed further below, are as follows:

  * •

_Fuzzy order relations_ [4] are generalizations of the conventional order
relations on the reals or the integer numbers: $\mathtt{SMALLER}$,
$\mathtt{EQUAL}$ and $\mathtt{GREATER}$. They enable a smooth transition
between these predicates and allow for expressing, for instance, that a number
$x$ is smaller than $y$ _to a certain degree_ , while to some degree these
numbers are also considered as being equal. Here, the $\mathtt{EQUAL}$
relation is understood as a kind of similarity relation that seeks to model
the “perceived equality” (instead of the strict mathematical equality).

  * •

_Scaling functions_ for modeling fuzzy equivalence relations [12]. For each
element $x$ of a linearly ordered domain $X$, a scaling function $s(\cdot)$
essentially expresses the degree $s(x)$ to which $x$ can be (or should be)
distinguished from its neighboring values. A measure of distance (or,
equivalently, of similarity) on $X$ can then be derived via accumulation of
local degrees of distinguishability.

  * •

_Fuzzy rank correlation_ [5, 17] generalizes conventional rank correlation on
the basis of fuzzy order relations, thereby combining properties of standard
rank correlation (such as Kendall’s tau) and numerical correlation measures
(such as Pearson correlation). Roughly, the idea is to penalize the inversion
of two items (later on called a _discordance_) depending on how dissimilar the
corresponding rank positions are: the more similar (less distinguishable) the
positions are according to the $\mathtt{EQUAL}$ relation, the smaller the
influence of the inversion on the rank correlation.

The rest of the paper is organized as follows. In the next two sections, we
briefly recall the basics of fuzzy order relations and fuzzy rank correlation,
respectively. Our weighted rank correlation measure, called _scaled gamma_ ,
is then introduced in Section 5, and related work is reviewed in Section 6. A
small experimental study is presented in Section 7, prior to concluding the
paper in Section 8.

##  2 Rank Correlation

Consider $N\geq 2$ paired observations
$\\{(x_{i},y_{i})\\}_{i=1}^{N}\subset\mathbb{X}\times\mathbb{Y}$ of two
variables $X$ and $Y$, where $\mathbb{X}$ and $\mathbb{Y}$ are two linearly
ordered domains; we denote

$\boldsymbol{x}=(x_{1},x_{2},\ldots,x_{N}),\enspace\boldsymbol{y}=(y_{1},y_{2},\ldots,y_{N})\enspace.$

In particular, the values $x_{i}$ (and $y_{i}$) can be real numbers
($\mathbb{X}=\mathbb{R}$) or rank positions
($\mathbb{X}=[N]\mathrel{\vbox{\hbox{\footnotesize.}\hbox{\footnotesize.}}}=\\{1,2,\ldots,N\\}$).
For example, $\boldsymbol{x}=(3,1,4,2)$ denotes a ranking of four items, in
which the first item is on position 3, the second on position 1, the third on
position 4 and the fourth on position 2.

The goal of a (rank) correlation measure is to capture the dependence between
the two variables in terms of their tendency to increase and decrease (their
position) in the same or the opposite direction. If an increase in $X$ tends
to come along with an increase in $Y$, then the (rank) correlation is
positive. The other way around, the correlation is negative if an increase in
$X$ tends to come along with a decrease in $Y$. If there is no dependency of
either kind, the correlation is (close to) 0.

###  2.1 Concordance and Discordance

Many rank correlation measures are defined in terms of the number $C$ of
_concordant_ , the number $D$ of _discordant_ , and the number $T$ of _tied_
data points. Let $\mathcal{P}=\\{(i,j)\,|\,1\leq i<j\leq N\\}$ denote the set
of ordered index pairs. We call a pair $(i,j)\in\mathcal{P}$ concordant,
discordant or tied depending on whether $(x_{i}-x_{j})(y_{i}-y_{j})$ is
positive, negative or 0, respectively. Thus, let us define the following
$N\times N$ relations:

$\displaystyle\mathcal{C}(i,j)$
$\displaystyle=\left\\{\begin{array}[]{cl}1&\text{ if
}\operatorname{sign}(\pi_{x}(i)-\pi_{x}(j))\operatorname{sign}(\pi_{y}(i)-\pi_{y}(j))=1\\\
0&\text{ otherwise}\end{array}\right.$ (3) $\displaystyle\mathcal{D}(i,j)$
$\displaystyle=\left\\{\begin{array}[]{cl}1&\text{ if
}\operatorname{sign}(\pi_{x}(i)-\pi_{x}(j))\operatorname{sign}(\pi_{y}(i)-\pi_{y}(j))=-1\\\
0&\text{ otherwise}\end{array}\right.$ (6) $\displaystyle\mathcal{T}_{x}(i,j)$
$\displaystyle=\left\\{\begin{array}[]{cl}1&\text{ if
}\operatorname{sign}(\pi_{x}(i)-\pi_{x}(j))=0,\,\operatorname{sign}(\pi_{y}(i)-\pi_{y}(j))\neq
0\\\ 0&\text{ otherwise}\end{array}\right.$ (9)
$\displaystyle\mathcal{T}_{y}(i,j)$
$\displaystyle=\left\\{\begin{array}[]{cl}1&\text{ if
}\operatorname{sign}(\pi_{x}(i)-\pi_{x}(j))\neq
0,\,\operatorname{sign}(\pi_{y}(i)-\pi_{y}(j))=0\\\ 0&\text{
otherwise}\end{array}\right.$ (12) $\displaystyle\mathcal{T}_{x,y}(i,j)$
$\displaystyle=\left\\{\begin{array}[]{cl}1&\text{ if
}\operatorname{sign}(\pi_{x}(i)-\pi_{x}(j))=0,\,\operatorname{sign}(\pi_{y}(i)-\pi_{y}(j))=0\\\
0&\text{ otherwise}\end{array}\right.$ (15) $\displaystyle\mathcal{T}(i,j)$
$\displaystyle=\left\\{\begin{array}[]{cl}1&\text{ if
}\operatorname{sign}(\pi_{x}(i)-\pi_{x}(j))\operatorname{sign}(\pi_{y}(i)-\pi_{y}(j))=0\\\
0&\text{ otherwise}\end{array}\right.$ (18)

Obviously,
$\mathcal{T}=\mathcal{T}_{x}\cup\mathcal{T}_{y}\cup\mathcal{T}_{x,y}$.

The number of concordant, discordant and tied pairs $(i,j)\in\mathcal{P}$ are
then obtained by summing the entries in the corresponding relations:

$\displaystyle C$
$\displaystyle=\sum_{(i,j)\in\mathcal{P}}\mathcal{C}(i,j)=\frac{1}{2}\sum_{i\in[N]}\sum_{j\in[N]}\mathcal{C}(i,j)$
$\displaystyle D$
$\displaystyle=\sum_{(i,j)\in\mathcal{P}}\mathcal{D}(i,j)=\frac{1}{2}\sum_{i\in[N]}\sum_{j\in[N]}\mathcal{D}(i,j)$
$\displaystyle T$
$\displaystyle=\sum_{(i,j)\in\mathcal{P}}\mathcal{T}(i,j)=\frac{1}{2}\sum_{i\in[N]}\sum_{j\in[N]}\mathcal{T}(i,j)-\frac{N}{2}$

Note that

$\mathcal{C}(i,j)+\mathcal{D}(i,j)+\mathcal{T}(i,j)=1$ (19)

for all $(i,j)\in\mathcal{P}$, and

$C+D+T=|\mathcal{P}|=\frac{N(N-1)}{2}\enspace.$ (20)

Well-known examples of rank correlation measures that can be expressed in
terms of the above quantities include Kendall’s tau [11]

$\tau=\frac{C-D}{N(N-1)/2}$ (21)

and Goodman and Kruskal’s gamma coefficient [9]

$\gamma=\frac{C-D}{C+D}\enspace.$ (22)

As will be detailed in the following sections, our basic strategy for
generalizing rank correlation measures such as $\gamma$ is to “fuzzify” the
concepts of concordance and discordance. Thanks to the use of fuzzy order
relations, we will be able to express that a pair $(i,j)$ is concordant or
discordant to a certain degree (between $0$ and $1$). Measures like (22) can
then be generalized in a straightforward way, namely by accumulating the
degrees of concordance and discordance, respectively, and putting them in
relation to each other.

###  2.2 Rank Correlation and Distance

Formally, there is a close connection between rank correlation and distance
measures on rankings. Indeed, many rank correlation measures are in fact
normalized versions of corresponding distance measures. For example,
Spearman’s rho is an affine transformation of the sum of squared rank
distances to the interval $[-1,+1]$, and Kendall’s tau is a similar
transformation of the Kendall distance, namely the sum of rank inversions. In
general, the relationship between a correlation $\operatorname{corr}$ and a
distance $d$ is expressed as follows:

$\operatorname{corr}(\pi_{x},\pi_{y})=1-\frac{2d(\pi_{x},\pi_{y})}{M}\,,$ (23)

where $M$ is a maximal distance. Desirable formal properties of a measure are
more naturally given in terms of distance. For example, the following axioms
were proposed by Kemmey and Snell [7]:

  * A1

The distance $d$ is a metric, i.e., the following properties hold for all
$\pi_{x},\pi_{y},\pi_{z}$:

    * –

$d(\pi_{x},\pi_{x})\geq 0$

    * –

$(\pi_{x}=\pi_{y})$ implies $d(\pi_{x},\pi_{y})=0$

    * –

$d(\pi_{x},\pi_{y})=0$ implies $\pi_{x}=\pi_{y}$

    * –

$d(\pi_{x},\pi_{y})=d(\pi_{x},\pi_{y})$

    * –

$d(\pi_{x},\pi_{z})\leq d(\pi_{x},\pi_{y})+d(\pi_{y},\pi_{z})$.

  * A2

Re-labeling: If $\pi_{x}^{\prime}$ results from $\pi_{x}$ by a permutation of
the items, and $\pi_{y}^{\prime}$ results from $\pi_{y}$ by the same
permutation, then $d(\pi_{x}^{\prime},\pi_{y}^{\prime})=d(\pi_{x},\pi_{y})$ .

  * A3

If two rankings $\pi_{x}$ and $\pi_{y}$ agree on all except a subset $S$ of
items, and these items form a segment in both rankings (i.e., they occupy a
sequence of consecutive positions), then $d(\pi_{x},\pi_{y})$ equals the
distance on the subset $S$ (i.e., the rankings $\pi_{x}$ and $\pi_{y}$
projected to $S$).

  * A4

The minimum positive distance is 1.

##  3 Fuzzy Relations

###  3.1 Fuzzy Equivalence

The notion of a fuzzy relation generalizes the standard notion of a
mathematical relation by allowing one to express “degrees of relatedness”.
Formally, a (binary) fuzzy relation on a set $\mathbb{X}$ is characterized by
a membership function

$\mathcal{E}:\mathbb{X}\times\mathbb{X}\longrightarrow\left[0,1\right]\,.$

For all $x,y\in\mathbb{X}$, $\mathcal{E}(x,y)$ is the degree to which $x$ is
related to $y$.

Recall that a conventional equivalence relation on a set $\mathbb{X}$ is a
binary relation that is reflexive, symmetric and transitive. For the case of a
fuzzy relation $\mathcal{E}$, these properties are generalized as follows:

  * •

reflexivity: $\mathcal{E}(x,x)=1$ for all $x\in\mathbb{X}$

  * •

symmetry: $\mathcal{E}(x,y)=\mathcal{E}(y,x)$ for all $x,y\in\mathbb{X}$

  * •

$\top$-transitivity:
$\top(\mathcal{E}(x,y),\mathcal{E}(y,z))\leq\mathcal{E}(x,z)$ for all
$x,y,z\in\mathbb{X}$

A fuzzy relation $\mathcal{E}$ having these properties is called a fuzzy
equivalence relation [5]. While the generalizations of reflexivity and
symmetry are rather straightforward, the generalization of transitivity
involves a triangular norm (t-norm) $\top$, which plays the role of a
generalized logical conjunction [13]. Formally, a function
$\top:\,[0,1]^{2}\longrightarrow[0,1]$ is a t-norm if it is associative,
commutative, monotone increasing in both arguments, and satisfies the boundary
conditions $\top(a,0)=0$ and $\top(a,1)=a$ for all $a\in[0,1]$. Examples of
commonly used t-norms include the minimum $\top(a,b)=\min(a,b)$ and the
product $\top(a,b)=ab$. To emphasize the role of the t-norm, a relation
$\mathcal{E}$ satisfying the above properties is also called a
$\top$-equivalence.

###  3.2 Fuzzy Ordering

The notion of an order relation $\leq$ is similar to that of an equivalence
relation, with the important difference that the former is antisymmetric while
the latter is symmetric. A common way to formalize antisymmetry is as follows:
$a\leq b$ and $b\leq a$ implies $a=b$. Note that this definition already
involves an equivalence relation, namely the equality $=$ of two elements.
Thus, as suggested by Bodenhofer [2], a fuzzy order relation can be defined on
the basis of a fuzzy equivalence relation. Formally, a fuzzy relation
$\mathcal{L}:\mathbb{X}\times\mathbb{X}\longrightarrow\left[0,1\right]$ is
called a _fuzzy ordering_ with respect to a t-norm $\top$ and a
$\top$-equivalence $\mathcal{E}$, for brevity $\top$-$\mathcal{E}$-$ordering$,
if it satisfies the following properties for all $x,y,z\in\mathbb{X}$:

  * •

$\mathcal{E}$-reflexivity: $\mathcal{E}(x,y)\leq\mathcal{L}(x,y)$

  * •

$\top$-$\mathcal{E}$-antisymmetry:
$\top(\mathcal{L}(x,y),\mathcal{L}(y,x))\leq\mathcal{E}(x,y)$

  * •

$\top$-transitivity:
$\top(\mathcal{L}(x,y),\mathcal{L}(y,z))\leq\mathcal{L}(x,z)$

Furthermore a $\top$-$\mathcal{E}$-ordering $\mathcal{L}$ is called strongly
complete if

$\max\big{(}\mathcal{L}(x,y),\mathcal{L}(y,x)\big{)}=1$

for all $x,y\in\mathbb{X}$. This is expressing that, for each pair of elements
$x$ and $y$, either $x\leq y$ or $y\leq x$ should be fully true.

A fuzzy relation $\mathcal{L}$ as defined above can be seen as a
generalization of the conventional “smaller or equal” on the real or the
integer numbers. What is often needed, too, is a “stricly smaller” relation
$<$. In agreement with the previous formalizations, a relation of that kind
can be defined as follows: A binary fuzzy relation $\mathcal{R}$ is called a
_strict fuzzy ordering_ with respect to a $\top$-norm and a $\top$-equivalence
$\mathcal{E}$, or strict $\top$-$\mathcal{E}$-$ordering$ for short, if it has
the following properties for all $x,x^{\prime},y,y^{\prime},z\in\mathbb{X}$
[5]:

  * •

irreflexivity: $\mathcal{R}(x,x)=0$

  * •

$\top$-transitivity:
$\top(\mathcal{R}(x,y),\mathcal{R}(y,z))\leq\mathcal{R}(x,z)$

  * •

$\mathcal{E}$-extensionality:
$\top(\mathcal{E}(x,x^{\prime}),\mathcal{E}(y,y^{\prime}),\mathcal{R}(x,y))\leq\mathcal{R}(x^{\prime},y^{\prime})$

###  3.3 Practical Construction

The above definitions provide generalizations $\mathcal{E}$, $\mathcal{L}$ and
$\mathcal{R}$ of the standard relations $=$, $\leq$ and $<$, respectively,
that exhibit reasonable properties and, moreover, are coherent with each
other. Practically, one may start by choosing an equivalence relation
$\mathcal{E}$ and a compatible t-norm $\top$, and then derive $\mathcal{L}$
and $\mathcal{R}$ from the corresponding $\top$-equivalence.

More specifically, suppose the set $\mathbb{X}$ to be a linearly ordered
domain, that is, to be equipped with a standard (non-fuzzy) order relation
$\leq$. Then, given a $\top$-equivalence $\mathcal{E}$ on $\mathbb{X}$, the
following relation is a coherent fuzzy order relation, namely a strongly
complete $\top$-$\mathcal{E}$-ordering:

$\mathcal{L}(x,y)=\left\\{\begin{array}[]{cl}1&\text{if }x\leq y\\\
\mathcal{E}(x,y)&\text{otherwise }\end{array}\right.$

Moreover, a strict fuzzy ordering $\mathcal{R}$ can be obtained from
$\mathcal{L}$ by

$\mathcal{R}(x,y)=1-\mathcal{L}(y,x)$ (24)

The relations thus defined have a number of convenient properties. In
particular, $\min(\mathcal{R}(x,y),\mathcal{R}(y,x))=0$ and

$\mathcal{R}(x,y)+\mathcal{E}(x,y)+\mathcal{R}(y,x)=1$ (25)

for all $x,y\in\mathbb{X}$. These properties can be interpreted as follows.
For each pair of elements $x$ and $y$, the unit mass splits into two parts: a
degree $a=\mathcal{E}(x,y)$ to which $x$ and $y$ are equal, and a degree $1-a$
to which either $x$ is smaller than $y$ or $y$ is smaller than $x$.

##  4 Fuzzy Relations on Rank Data

Since we are interested in generalizing rank correlation measures, the
underlying domain $\mathbb{X}$ is given by a set of rank positions $[N]$
(equipped with the standard $<$ relation) in our case. As mentioned before,
this domain could be equipped with fuzzy relations $\mathcal{E}$,
$\mathcal{L}$ and $\mathcal{R}$ by defining $\mathcal{E}$ first and deriving
$\mathcal{L}$ and $\mathcal{R}$ afterward. Note, however, that the number of
degrees of freedom in the specification of $\mathcal{E}$ is of the order
$O(N^{2})$, despite the constraints this relation has to meet.

###  4.1 Scaling Functions on Rank Positions

In order to define fuzzy relations even more conveniently, while emphasizing
the idea of weighing the importance of rank positions at the same time, we
leverage the concept of a _scaling function_ as proposed by Klawonn [12].
Roughly speaking, a scaling function
$w:\,\mathbb{X}\longrightarrow\mathbb{R}_{+}$ specifies the dissimilarity of
an element $x$ from its direct neighbor elements, and the dissimilarity
between any two elements $x$ and $y$ is then obtained via integration of the
local dissimilarities along the chain from $x$ to $y$. In our case, a scaling
function can be defined as a mapping $w:\,[N-1]\longrightarrow[0,1]$ or,
equivalently, as a vector

$\vec{w}=\Big{(}w(1),w(2),\ldots,w(N-1)\Big{)}\in[0,1]^{N-1}\enspace.$ (26)

Here, $w(n)$ can be interpreted as the degree to which the rank positions $n$
and $n+1$ are distinguished from each other; correspondingly, $1-w(n)$ can be
seen as the degree to which these two positions are considered to be equal.
From the local degrees of distinguishability, a global distance function is
derived on $\mathbb{X}$ by defining

$d(x,y)=\min\left(1,\sum_{i=\min(x,y)}^{\max(x,y)-1}w(i)\right)\enspace.$ (27)

Put in words, the distance between $x$ and $y$ is the sum of the degrees of
distinguishability between them, thresholded at the maximal distance of 1. In
principle, accumulations of the degrees of distinguishability other than the
sum are of course conceivable. For example, the maximum could be used as well:

$d(x,y)=\max\Big{\\{}w(i)\,|\,i\in\\{\min(x,y),\ldots,\max(x,y)-1\\}\Big{\\}}\enspace.$
(28)

In general, $d(x,y)$ is supposed to define a pseudo-metric on $\mathbb{X}$.
Under this condition, it can be shown that the fuzzy relation $\mathcal{E}$
defined as

$\mathcal{E}(x,y)=1-d(x,y)$

for all $x,y\in\mathbb{X}$ is a $\top_{L}$-equivalence, where $\top_{L}$ is
the Łukasiewicz $t$-norm $\top_{L}(a,b)=\max\\{0,a+b-1\\}$ [3]. Relations
$\mathcal{L}$ and $\mathcal{R}$ can then be derived from $\mathcal{E}$ as
described in Section 3.3. In particular, we obtain

$\mathcal{R}(x,y)=\left\\{\begin{array}[]{cl}d(x,y)&\text{ if }x<y\\\ 0&\text{
otherwise }\end{array}\right.$

According to our discussion so far, the only remaining degree of freedom is
the scaling function $s$. Obviously, this function can also be interpreted as
a _weighing function_ : the more distinguishable a position $n$ from its
neighbor positions, i.e., the larger $w(n-1)$ and $w(n)$, the higher the
importance of that position.

Figure 1: Example of a scaling function.

An example of a scaling function for $N=12$ is shown in Figure 1. This
function puts more emphasis on the top and the bottom ranks and less on the
middle part. According to (27), the distinguishability between the positions
$4$ and $7$ is $d(4,7)=0.4+0.2+0.2=0.8$ (sum of the weights $w(i)$ in the
shaded region). Thus, $4$ is strictly smaller than $7$ to the degree of
$\mathcal{R}(4,7)=0.8$, while both positions are considered equal to the
degree $\mathcal{E}(4,7)=0.2$.

Note that, with $w(i)=\llbracket i<k\rrbracket$, we also cover the top-$k$
scenario as a special case. Here, the standard $<$ relation is recovered for
all elements on the first $k$ positions, whereas the remaining positions are
considered as fully equivalent, i.e., these elements form an equivalence class
in the standard sense.

##  5 Weighted Rank Correlation

Our approach to generalizing rank correlation measures is based on the
“fuzzificiation” of the relations (3–18) and, correspondingly, the number of
concordant, discordant and tied item pairs. The tools that are needed to do so
have already been introduced in the previous sections. In particular, suppose
a fuzzy equivalence relation $\mathcal{E}$ and a “strictly smaller” relation
$\mathcal{R}$ to be derived from a scaling function $w$ on $\mathbb{X}$, based
on the procedure outlined above. For notational convenience, we assume the
same scaling function (and hence the same relations) to be used on both
domains $\mathbb{X}$ and $\mathbb{Y}$. In principle, however, different
functions $w_{X}$ and $w_{Y}$ (and hence relations $\mathcal{E}_{X}$,
$\mathcal{R}_{X}$ and $\mathcal{E}_{Y}$, $\mathcal{R}_{Y}$) could be used.

Now, according to (3), a pair $(i,j)\in\mathcal{P}$ is concordant if both
$x_{i}$ is (strictly) smaller than $x_{j}$ and $y_{i}$ is smaller than
$y_{j}$, or if $x_{j}$ is smaller than $x_{i}$ and $y_{j}$ is smaller than
$y_{i}$. Using our fuzzy relation $\mathcal{R}$ and a t-norm $\top$ as a
generalized conjunction, this can be expressed as follows:

$\tilde{\mathcal{C}}(i,j)=\top\big{(}\mathcal{R}(x_{i},x_{j}),\mathcal{R}(y_{i},y_{j})\big{)}+\top\big{(}\mathcal{R}(x_{j},x_{i}),\mathcal{R}(y_{j},y_{i})\big{)}$
(29)

The discordance relation can be expressed analogously:

$\tilde{\mathcal{D}}(i,j)=\top\big{(}\mathcal{R}(x_{i},x_{j}),\mathcal{R}(y_{j},y_{i})\big{)}+\top\big{(}\mathcal{R}(x_{j},x_{i}),\mathcal{R}(y_{i},y_{j})\big{)}$
(30)

Finally, the degree to which $(i,j)$ is tied is given by

$\tilde{\mathcal{T}}(i,j)=\bot\big{(}\mathcal{E}(x_{i},x_{j}),\mathcal{E}(y_{i},y_{j})\big{)}\enspace,$

where $\bot$ is the t-conorm associated with $\top$ (i.e.,
$\bot(u,v)=1-\top(1-u,1-v)$), serving as a generalized logical disjunction.
Generalizing (19), the three degrees sum up to 1, i.e.,

$\tilde{\mathcal{C}}(i,j)+\tilde{\mathcal{D}}(i,j)+\tilde{\mathcal{T}}(i,j)\equiv
1\enspace,$ (31)

and either $\tilde{\mathcal{C}}(i,j)=0$ or $\tilde{\mathcal{D}}(i,j)=0$. In
other words, a pair $(i,j)$ that has originally been concordant (discordant)
will remain concordant (discordant), at least to some extent. However, since
$\mathcal{E}$ may introduce a certain indistinguishability between the
positions $x_{i}$ and $x_{j}$ or the positions $y_{i}$ and $y_{j}$, the pair
could also be considered as a partial tie.

Given the above fuzzy relations, the number of concordant, discordant and tied
data points can be obtained as before, namely by summing over all ordered
pairs $(i,j)\in\mathcal{P}$:

$\tilde{C}=\sum_{(i,j)\in\mathcal{P}}\tilde{\mathcal{C}}(i,j)\enspace,\quad\tilde{D}=\sum_{(i,j)\in\mathcal{P}}\tilde{\mathcal{D}}(i,j)\enspace,\quad\tilde{T}=\sum_{(i,j)\in\mathcal{P}}\tilde{\mathcal{T}}(i,j)\enspace.$

According to (31),

$\tilde{C}+\tilde{D}+\tilde{T}=|\mathcal{P}|=\frac{N(N-1)}{2}\enspace,$

which generalizes (20). Using these quantities, rank correlation measures
expressed in terms of the number of concordant and discordant pairs can be
generalized in a straightforward way. In particular, a generalization of the
gamma coefficient (22) is obtained as

$\tilde{\gamma}=\frac{\tilde{C}-\tilde{D}}{\tilde{C}+\tilde{D}}\enspace.$ (32)

It is worth mentioning that the weighted rank correlation measure thus defined
exhibits a number of desirable formal properties, which it essentially
inherits from the general fuzzy extension of the gamma coefficient; we refer
to [17], in which these properties are analyzed in detail.

##  6 Dealing with Ties

The basic principle of our approach to positional weighing of rank correlation
measures is the idea of discounting concordances and discordances of pairs of
items by turning them into _partial ties_ : If the positions occupied by the
items are considered as (partly) equivalent according to the underlying
equivalence relation on ranks, a concordance (discordance) is turned into a
partial concordance (discordance) and a partial tie. The way in which ties are
dealt with by the subsequently applied rank correlation measure is hence very
critical:

  * •

To achieve the desired weighting effect, the effect of ties on the overall
correlation should be less strong than the effect of concordances or
discordances.

  * •

In this regard, it is important to realize that ties are indeed treated
differently by different measures, which is also due to different
interpretations of ties.

  * •

This issue becomes especially critical if we seek to apply our approach to
rankings that may already contain ties (which we excluded so far). In this
case, there are two types of ties: “real ties” contained in the original
rankings and “artificial ties” resulting from our transformation. These two
types should arguably be distinguished.

###  6.1 Notation

As before, we consider data in the form of objects
$\boldsymbol{x}=(x_{1},x_{2},\ldots,x_{N})\in\mathbb{X}=[N]\mathrel{\vbox{\hbox{\footnotesize.}\hbox{\footnotesize.}}}=\\{1,2,\ldots,N\\}$.
Now, however, we do not necessarily assume that all entries $x_{i}$ are
distinct. Instead, if $x_{i}=x_{j}$ for $i\neq j$, we consider the two items
to be tied. For example, denoting the items to be ranked by $A,B,C$, etc.,
$\vec{x}=(2,1,3,2,1)$ represents the ranking $\\{B,E\\}\succ\\{A,D\\}\succ C$
for items $\\{A,B,C,D,E\\}$, which we can also write as $BE|AD|C$.

As before, the preference relation $\succ$ on items is derived from the
standard $<$ relation on $\mathbb{X}$. The set of items with the same rank
form a tie group or equivalence class. We denote by $[n]=\\{i\,|\,x_{i}=n\\}$
the set of items (or rather their indices) put on position $n$. The preference
$\succ$ thus defined is not only irreflexive, asymmetric, and transitive, but
also satisfies transitivity of incomparability $\sim$, where $A\sim B$ if
neither $A\succ B$ nor $B\succ A$: If $A\sim B$ and $B\sim C$, then $A\sim C$.
In the literature, such relations are called strict weak order, or sometimes
“bucket orders” (with the idea that all items within an equivalence class are
put in a bucket).

###  6.2 Interpretation of Ties

There are two main interpretations of ties:

  * •

Indifference ($\sim_{E}$): First, two items can be tied because they are
indeed considered as equal ($\sim$ is interpreted as $\sim_{E}$). For example,
a decision maker might be indifferent between two choice alternatives, because
she likes them both the same.

  * •

Incomparability ($\sim_{U}$): Second, two items could be tied in a ranking
because their true relationship is unknown, i.e., they are considered as
incomparable due to a lack of information ($\sim$ is interpreted as
$\sim_{U}$).

It is important to notice that, according to the first interpretation, a tie
$A\sim B$ directly refers to properties of the underlying items $A$ and $B$
(it has an _ontic_ meaning), whereas the second interpretation is of
_epistemic_ nature and represents _knowledge_ about the relation between the
items: $A\sim B$ does not mean that $A$ and $B$ are equally preferred.
Instead, given the information at hand, we simply cannot compare them, i.e.,
we cannot determine the “true” relationship between $A$ and $B$. Also note
that, regarding this “ground truth”, different assumptions could be made: We
may assume that the ground truth is a total order, i.e., either $A\succ B$ or
$B\succ A$, we only do not know (yet). Or, we may even allow that the ground
truth is an indifference between $A$ and $B$. In this case, the two different
interpretations of ties would co-exist: $A\sim_{U}B$ means that $A\succ B$ or
$B\succ A$ or $A\sim_{E}B$. Obviously, this co-existence will lead to
additional complications. Therefore, whenever the interpretation of ties in
terms of incomparability is adopted, we assume that the ground truth is a
total order (ranking without indifference).

###  6.3 Treatment of Ties

Naturally, the different interpretations of ties call for different ways of
treating them in the context of rank correlation. For example, $A\sim_{E}B$ in
two rankings suggests a strong coherence, because both rankings agree on the
equality of $A$ and $B$, and this coherence should be rewarded by a
correlation measure. In the case of incomparability, the situation is slightly
different: $A\sim_{U}B$ in two observed rankings could indicate both a
concordance and a discordance in the underlying ground truth rankings. In this
situation, it might be natural to express _uncertainty_ about the actual
correlation.

More generally, the interpretation of ties in terms of indifference calls for
an extended scoring scheme, which, in addition to (scores for) concordance and
discordance also considers comparisons involving ties:

|  $\pi_{y}(i)<\pi_{y}(j)$ |  $\pi_{y}(i)=\pi_{y}(j)$ |  $\pi_{y}(i)>\pi_{y}(j)$
---|---|---|---
$\pi_{x}(i)<\pi_{x}(j)$ | $s_{1,1}$ | $s_{1,2}$ | $s_{1,3}$
$\pi_{x}(i)=\pi_{x}(j)$ | $s_{2,1}$ | $s_{2,2}$ | $s_{2,3}$
$\pi_{x}(i)>\pi_{x}(j)$ | $s_{3,1}$ | $s_{3,2}$ | $s_{3,3}$

Naturally, the highest scores (rewards) should be found on the diagonal,
whereas $s_{1,3}$ and $s_{3,1}$ should be lowest; the values
$s_{1,2},s_{2,1},s_{2,3},s_{3,2}$ should be in-between, because they can be
seen as half a concordance and half a discordance.

In the case of incomparability, it is natural to compare two rankings $\pi$
and $\pi^{\prime}$ by considering the set of possible values of correlation
between the underlying (but unknown) ground truth rankings, that is, the set
of extensions $E(\pi)$ and $E(\pi^{\prime})$. For a correlation measure
$\operatorname{corr}$, which only needs to be defined on complete rankings,
this set if given by

$\Big{\\{}\operatorname{corr}(\tau,\tau^{\prime})\,|\,\tau\in
E(\pi),\,\tau^{\prime}\in E(\pi^{\prime})\Big{\\}}\,.$

To obtain a single value, the values in this set could be aggregated in one
way or the other, for example in terms of the minimum (reflecting a
pessimistic attitude), the maximum (reflecting an optimistic attitude), or an
average.

##  7 Related Work

Weighted versions of rank correlation measures have not only been studied in
statistics but also in other fields, notably in information retrival [20, 10,
6, 15]. Most of them are motivated by the idea of giving a higher weight to
the top-ranks: in information retrieval, important documents are supposed to
appear in the top, and a swap of important documents should incur a higher
penalty than a swap of unimportant ones.

Kaye [10] introduced a weighted, non-symmetric version of Spearman’s rho
coefficient. Costa and Soares [6] proposed a symmetric weighted version of
Spearman’s coefficient resembling the one of Kaye. Another approach, based on
average precision and called AP correlation, was introduced by Yilmaz et al.
[20]. Maturi and Abdelfattah [15] define weighted scores $W_{i}=w^{i}$ with
$w\in\left(0,1\right)$ and compute the Pearson correlation coefficient on
these scores. All four measures give higher weight to the top ranks.

Two more flexible measures, not restricted to monotone decreasing weights,
have been proposed by Shieh [18] and Kumar and Vassilivitskii [14]. In the
approach of Shieh [18], a weight is manually given to every occurring
concordance or discordance through a symmetric weight function
$w:[N]\times[N]\longrightarrow\mathbb{R}_{+}$:

$\tau_{w}=\frac{\sum_{i<j}{w_{ij}C_{ij}}-\sum_{i<j}{w_{ij}D_{ij}}}{\sum_{i<j}{w_{ij}}}=\frac{\sum_{i<j}{w_{ij}(C_{ij}}-D_{ij})}{\sum_{i<j}{w_{ij}}}.$
(33)

The input parameter for $w$ are the ranks of a reference ranking $\pi_{ref}$,
which is assumed to be the natural order $\left(1,2,3,\ldots,N\right)$.
Therefore, this approach is not symmetric. To handle the quadratic number of
weights, Shieh proposed to define them as $w_{ij}=v_{i}v_{j}$ with $v_{i}$ the
weight of rank $i$.

Kumar and Vassilivitskii [14] introduce a generalized version of Kendall’s
distance. Originally, they proposed three different weights: element weights,
position weights, and element similarities. The three weights are defined
independently of each other, and each of them can be used by its own for
weighting discordant pairs. Here, we focus on the use of position weights.
Like in our approach, Kumar and Vassilivitskii define $N-1$ weights
$\delta_{i}\geq 0$, which are considered as costs for swapping two elements on
adjacent positions $i+1$ and $i$. The accumulated cost of changing from
position $1$ to $i\in\\{2,\ldots,N\\}$ is $p_{i}=\sum_{j=1}^{i-1}\delta_{j}$,
with $p_{1}=0$. Moreover,

$\bar{p}_{i}(\pi_{1},\pi_{2})=\frac{p_{\pi_{1}(i)}-p_{\pi_{2}(i)}}{\pi_{1}(i)-\pi_{2}(i)}$
(34)

is the average cost of moving element $i$ from position $\pi_{1}(i)$ to
position $\pi_{2}(i)$; if $\pi_{1}(i)=\pi_{2}(i)$ then $\bar{p}_{i}=1$.The
weighted discordance of a pair $(i,j)$ is then defined in terms of the product
of the average costs for index $i$ and $j$:

$\hat{D}_{\delta}(i,j)=\begin{cases}\bar{p}_{i}(\pi_{1},\pi_{2})\bar{p}_{j}(\pi_{1},\pi_{2})&\text{if}\
(i,j)\ \text{is discordant}\\\ 0&\text{otherwise}\end{cases}\ .$ (35)

Finally, the weighted Kendall distance $K_{\delta}$ is given by

$K_{\delta}=\tilde{D}_{\delta}\
=\sum_{i=1}^{N-1}{\sum_{i+1}^{N}{\hat{D}_{\delta}(i,j)}}\ .$ (36)

Note that (36) is indeed a distance and not a correlation measure. To ease
comparison with $\tau_{\omega}$ and $\tilde{\gamma}$, we can define

$\hat{C}_{\delta}(i,j)=\begin{cases}\bar{p}_{i}(\pi_{1},\pi_{2})\bar{p}_{j}(\pi_{1},\pi_{2})&\text{if}\
(i,j)\ \text{is concordant}\\\ 0&\text{otherwise}\end{cases}$ (37)

as the weighted concordance of a pair $(i,j)$, and finally another weighted
version of gamma:

$\tilde{\gamma}_{\delta}=\frac{\tilde{C}_{\delta}-\tilde{D}_{\delta}}{\tilde{C}_{\delta}+\tilde{D}_{\delta}}\enspace.$

Specific emphasis on the handling of ties is put in the recent paper [16]. The
authors study the problem of rank aggregation, i.e., aggregating preferences
into a consensus ranking, taking into account that swapping elements on the
top might be more important than swapping elements at the bottom. To this end,
they propose a position weighted rank correlation coefficient for rankings
with ties, which combines the structure of Emond and Mason’s measure [7] with
the class of weighted Kemeny-Snell distances.

##  8 Summary and Conclusion

We introduced a new approach to weighted rank correlation based on fuzzy order
relations, as well as a concrete measure called _scaled gamma_. The latter
allows for specifying the importance of rank positions in a quite flexible and
convenient way by means of a scaling function. Thanks to the underlying formal
foundation, such a scaling function immediately translates into a concrete
version of our measure, in which the rank positions are processed within an
appropriate weighting scheme.

In future work, we plan to analyze the usefulness of our new approach to rank
correlation by means of experimental studies, both with synthetic and real
data. More specifically, we plan to compare our measure with others on various
applications and data analysis tasks, such as ranking or clustering. Moreover,
let us again highlight that our extension of gamma is actually not a single
measure but a family of measures, which is parameterized by the weight
function $w$ as well as the generalized logical conjunction (t-norm) used to
define concordance and discordance. While the former will typically be
specified as an external parameter by the user, the (fuzzy) logical operators
offer an interesting degree of freedom that could be used to optimally adapt
the measure to the application at hand. Again, this is an interesting
direction for future work. Finally, going beyond the gamma coefficient, we
also intend to apply our generalization to other rank correlation measures.

A Python implementation111This implementation was provided by Hoang Cong
Thanh, whose contribution is kindly acknowledged by the authors. of our
weighted rank correlation measure can be found here:
https://github.com/KIuML/weighted_rank_correlation

## References

  * [1] R. Balasubramaniyan, E. Hüllermeier, N. Weskamp, and J. Kämper.  Clustering of gene expression data using a local shape-based similarity measure.  Bioinformatics, 21(7):1069–1077, 2005.
  * [2] U. Bodenhofer.  A similarity-based generalization of fuzzy orderings preserving the classical axioms.  Int. J. Uncertainty, Fuzziness and Knowledge-based Systems, 8(5):593–610, 2000.
  * [3] U. Bodenhofer.  Representations and constructions of similarity-based fuzzy orderings.  Fuzzy Sets and Systems, 137:113–136, 2003.
  * [4] U. Bodenhofer and M. Demirci.  Strict fuzzy orderings with a given context of similarity.  Int. J. of Uncertainty, Fuzziness and Knowledge-Based Systems, 16(2):147–178, 2008.
  * [5] U. Bodenhofer and F. Klawonn.  Robust rank correlation coefficients on the basis of fuzzy orderings: Initial steps.  Mathware & Soft Computing, 15:5–20, 2008.
  * [6] J. Pinto da Costa and C. Soares.  A weighted rank measure of correlation.  Australian & New Zealand Journal of Statistics, 47(4):515–529, 2005.
  * [7] E.J. Emond and D.W. Mason.  A new rank correlation coefficient with application to the consensus ranking problem.  Journal of Multi-Criteria Decision Analysis, 11:17–28, 2002.
  * [8] R. Fagin, R. Kumar, and D. Sivakumar.  Comparing top k lists.  SIAM Journal on Discrete Mathematics, 17(1):134–160, 2003.
  * [9] L.A. Goodman and W.H. Kruskal.  Measures of Association for Cross Classifications.  Springer-Verlag, New York, 1979.
  * [10] D. Kaye.  A weighted rank correlation coefficient for the comparison of relevance judgements.  Journal of Documentation, 29(4):380–389, 1973.
  * [11] M.G. Kendall.  Rank correlation methods.  Charles Griffin, London, 1955.
  * [12] F. Klawonn.  Fuzzy sets and vague environment.  Fuzzy Sets and Systems, 66:207–221, 1994.
  * [13] E.P. Klement, R. Mesiar, and E. Pap.  Triangular Norms.  Kluwer Academic Publishers, 2002.
  * [14] R. Kumar and S. Vassilvitskii.  Generalized distances between rankings.  In Proc. WWW, 19. International Conference on World Wide Web, pages 571–580, 2010.
  * [15] T.A. Maturi and E.H. Abdelfattah.  A new weighted rank correlation.  Journal of Mathematics and Statistics, 4(4):226, 2008.
  * [16] A. Plaia, S. Buscemi, and M. Sciandra.  Consensus among preference rankings: a new weighted correlation coefficient for linear and weak orderings.  Advances in Data Analysis and Classification, 15:1015–1037, 2021\.
  * [17] M. Dolorez Ruiz and E. Hüllermeier.  A formal and empirical analysis of the fuzzy gamma rank correlation coefficient.  Information Sciences, 206:1–17, 2012.
  * [18] G.S. Shieh.  A weighted Kendall’s tau statistic.  Statistics & Probability Letters, 39(1):17–24, 1998.
  * [19] C. Spearman.  The proof and measurement for association between two things.  Amer. Journal of Psychology, 15:72–101, 1904.
  * [20] E. Yilmaz, J.A. Aslam, and S. Robertson.  A new rank correlation coefficient for information retrieval.  In Proc. 31st Annual International ACM SIGIR Conference on Research and Development in Information Retrieval, pages 587–594. ACM, 2008\.