# Dynamic Ranking with the BTL Model: A Nearest Neighbor based Rank Centrality Method

Eglantine Karlé, Hemant Tyagi

<figure><img src="image_1.png" /></figure>

## To cite this version:

Eglantine Karlé, Hemant Tyagi. Dynamic Ranking with the BTL Model: A Nearest Neighbor based Rank Centrality Method. Journal of Machine Learning Research, 2023, 24 (269), pp.1-57. hal-03519271

HAL Id: hal-03519271

https://inria.hal.science/hal-03519271

Submitted on 10 Jan 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>
# Dynamic Ranking with the BTL Model: A Nearest Neighbor based Rank Centrality Method

Eglantine Karlé *¹ and Hemant Tyagi †¹

¹Inria, Univ. Lille, CNRS, UMR 8524 - Laboratoire Paul Painlevé, F-59000.

September 29, 2021

## Abstract

Many applications such as recommendation systems or sports tournaments involve pairwise comparisons within a collection of *n* items, the goal being to aggregate the binary outcomes of the comparisons in order to recover the latent strength and/or global ranking of the items. In recent years, this problem has received significant interest from a theoretical perspective with a number of methods being proposed, along with associated statistical guarantees under the assumption of a suitable generative model.

While these results typically collect the pairwise comparisons as one comparison graph *G*, however in many applications – such as the outcomes of soccer matches during a tournament – the nature of pairwise outcomes can evolve with time. Theoretical results for such a dynamic setting are relatively limited compared to the aforementioned static setting. We study in this paper an extension of the classic BTL (Bradley-Terry-Luce) model for the static setting to our dynamic setup under the assumption that the probabilities of the pairwise outcomes evolve smoothly over the time domain [0, 1]. Given a sequence of comparison graphs $(G_{t'})_{t' \in \mathcal{T}}$ on a regular grid $\mathcal{T} \subset [0, 1]$, we aim at recovering the latent strengths of the items $w_t \in \mathbb{R}^n$ at any time $t \in [0, 1]$. To this end, we adapt the Rank Centrality method – a popular spectral approach for ranking in the static case – by locally averaging the available data on a suitable neighborhood of *t*. When $(G_{t'})_{t' \in \mathcal{T}}$ is a sequence of Erdös-Renyi graphs, we provide non-asymptotic $\ell_2$ and $\ell_\infty$ error bounds for estimating $w_t^*$ which in particular establishes the consistency of this method in terms of *n*, and the grid size $|\mathcal{T}|$. We also complement our theoretical analysis with experiments on real and synthetic data.

## 1 Introduction

Ranking problems arise in many domains such as sports tournaments [5, 2], recommendation systems [12] and even in the study of systems in biology [4, 28]. Typically, the data associated with these problems are binary outcomes of pairwise comparisons. For example in a sports tournament, the data comprises of the result of each game played between two teams or players (team *i* beat team *j*). The aim then is to recover a ranking of the items from the pairwise outcomes, often by estimating a score for each item, that represents its latent strength/quality. We remark that comparison data are frequently preferred to the attribution of a score by the users due to the sensitivity

<footnote>*Email: eglantine.karle@inria.fr</footnote>

<footnote>†Email: hemant.tyagi@inria.fr</footnote>
to bias of the latter. Indeed, the evaluation of the items and the strategy of scoring usually varies considerably between users and thus lead to rather inconsistent scores [12].

The collection of pairs of items which are compared can be naturally represented as an undirected graph $G = ([n], \mathcal{E})$ where the $n$ vertices are the items to rank, and an edge $\{i, j\} \in \mathcal{E}$ represents a comparison between items $i$ and $j$. Note that in order to recover the underlying ranking, the graph $G$ needs necessarily needs to be connected. Indeed, if this is not the case, then there will exist two connected components of $G$ such that there is no information available to rank items in one component with respect to the other.

**Classic BTL model.** A classic statistical model for ranking problems was introduced by Bradley and Terry [3] and completed by Luce [17]; this is now referred to as the BTL (Bradley-Terry-Luce) model. Denoting $w_i^* > 0$ to be the strength of item $i$, it posits that the probability that item $j$ beats item $i$, namely $y_{ij}^*$, is given by $\frac{w_j^*}{w_i^*+w_j^*}$. The data is then available in the form of independent Bernoulli random variables paramaterized by $y_{ij}^*$ for each $\{i, j\} \in \mathcal{E}$. This model has been studied extensively through many applications, such as sports tournaments [5], measurements of pain [18] or the estimation of the risk of car-crashes [14]. Moreover, numerous methods have been proposed to analyse the model theoretically with the goal of estimating the latent strengths. Two methods which have received particular attention are the maximum likelihood estimator (MLE) [3, 23], and a spectral method called Rank Centrality introduced by Neghaban et al. [23]. Both these methods have been analyzed in detail over the years for general graph topologies (eg., [23]) and also for the special case of Erdös-Renyi graphs $G \sim \mathcal{G}(n, p)$ (eg., [7, 6]). For instance, let us denote $\pi^* = \frac{w^*}{\|w^*\|_1}$ to be the normalized vector of true scores, and $\hat{\pi}$ its estimate returned by Rank Centrality. Then assuming that each comparison between two items has been performed $L$ times (for $L$ large enough) and $p \gtrsim \frac{\log n}{n}$, it holds with high probability (w.h.p) that [7, 6]

$$ \frac{\|\hat{\pi} - \pi^*\|_2}{\|\pi^*\|_2} \lesssim \frac{1}{\sqrt{npL}} \quad \text{and} \quad \frac{\|\hat{\pi} - \pi^*\|_\infty}{\|\pi^*\|_\infty} \lesssim \sqrt{\frac{\log n}{npL}}. $$

Here $\lesssim, \gtrsim$ hide absolute constants, see notations in Section 2. Analogous bounds hold for the MLE as well, and moreover, these bounds are optimal [7, 6]. The $\ell_\infty$ bounds are particularly useful as they readily lead to conditions for exact recovery of the true ranks of the items. A weaker requirement is the identification of the set of top-$K$ ranked items for some small $K$, see [7].

**Dynamic ranking.** The aforementioned setup considers the 'static' case where only one comparison graph is available. However, in many applications, the latent strengths of the items evolves over time. Indeed, a user preference changes with time in recommendation systems, as does the level of a player or a team in sports tournaments. This dynamic setting has received significantly less attention than its static counterpart and has mostly been studied with a focus on applications, such as sports tournaments [5, 11, 21]. It has rarely been analysed theoretically in the past however, although some results exist for a state-space generalization of the BTL model [9, 10, 19] and for a Bayesian framework [10, 16].

**Dynamic BTL.** The focus of this paper is to consider a dynamic version of the classical BTL model for ranking. In this *dynamic* scenario, we are given pairwise comparisons on a discrete set of time instants $T$, hence leading to a sequence of comparison graphs $(G_{t'})_{t' \in T}$. Moreover, denoting $w_{t,i}^* > 0$ to be the strength of item $i$ at time $t$, the probability that item $j$ beats item $i$ at time
$t' \in \mathcal{T}$ is now given by $\frac{w_{t',j}^*}{w_{t',i}^*+w_{t',j}^*}$. This model is described in detail in Section 2.1. The goal is to recover the underlying ranking at any given time instant $t$, where $t$ does not necessarily belong to $\mathcal{T}$. Furthermore, we make the following remarks.

* Suppose $t$ lies in $\mathcal{T}$, then we could consider only using the pairwise comparison data corresponding to $G_t$, along with a suitable algorithm for the static case. However, $G_t$ is not necessarily going to be connected; indeed, all the comparison graphs can now be disconnected and be very sparse. Hence this approach will fail in such scenarios. This suggests that it is in some sense crucial that we utilize the pairwise comparison information across different time instants in $\mathcal{T}$ in order to get meaningful rank recovery at time $t$. This approach can of course be relevant even if $G_{t'}$ was connected for all $t' \in \mathcal{T}$.

* Very recently, a couple of theoretical results were obtained for the dynamic ranking problem. Li and Wakin [15] include in their model the observed difference of scores (for example the score of a football game) which clearly contains richer information in comparison to simple binary outcomes. Bong et al. [2] considered an adaptation of the BTL model to the dynamic setting similar to what has been described above. They first combine pairwise comparison information by applying kernel smoothing to the entire dataset, and then use the MLE for the smoothed data to recover the ranks. We discuss these results in more detail in Section 7.

**Our contributions.** We provide in this paper the following contributions.

* We consider a model adapted to the dynamic ranking setting which involves using the BTL model at each time point, but assumes that the pairwise outcome probabilities are Lipschitz functions of time. In this setting, we propose an adaptation of the Rank Centrality method, namely Dynamic Rank Centrality (Algorithm 1) to recover the ranks at any given time $t$ by averaging the data belonging to a suitable neighborhood of $t$.

* We provide a detailed theoretical analysis of the previous algorithm by deriving $\ell_2$ and $\ell_\infty$ error bounds for estimating the latent strength vector $w_t^* = (w_{t,1}^*, \dots, w_{t,n}^*)^\top$ at any given time instant $t \in [0, 1]$. In particular, assuming $\mathcal{T} \subset [0, 1]$ is a uniform grid ($|\mathcal{T}| = T + 1$) and the comparison graphs $(G_{t'})_{t' \in \mathcal{T}}$ are Erdös-Renyi, we show (see Corollaries 1,2) for $T$ large enough that the error rates are of the order $O(T^{-1/3})$ w.h.p. This result holds in very sparse regimes where each $G_{t'}$ could be potentially disconnected. Moreover, this rate matches the classical pointwise estimation error rate for Lipschitz functions.

* Finally, we perform extensive experiments on synthetic data which validates our theoretical findings. We also evaluate our method on a real data set and show that it performs well.

**Outline of the paper.** After presenting the dynamic BTL model and our algorithm in Section 2, we summarize our main theoretical results for the $\ell_2$ and $\ell_\infty$ estimation errors in Section 3. We then describe in detail the $\ell_2$ analysis in Section 4, and the $\ell_\infty$ analysis in Section 5. Finally, we conduct experiments on both synthetic and real data in Section 6, and discuss other closely related work in Section 7.

# 2 Problem setup and algorithm

**Notation** For any probability vector $\pi \in \mathbb{R}^n$ with strictly positive entries, we define the vector norm $\|x\|_\pi = \sqrt{\sum_{i=1}^n \pi_i x_i^2}$. For a matrix $A$, the corresponding induced matrix norm is then defined
as $\|A\|_{\pi} = \sup_{\|x\|_{\pi}=1} \|x^{\top}A\|_{\pi}$. Note that some simple inequalities follow from these definitions.

$$ \sqrt{\pi_{\min}} \|x\|_2 \le \|x\|_{\pi} \le \sqrt{\pi_{\max}} \|x\|_2 \quad \text{and} \quad \sqrt{\frac{\pi_{\min}}{\pi_{\max}}} \|A\|_2 \le \|A\|_{\pi} \le \sqrt{\frac{\pi_{\max}}{\pi_{\min}}} \|A\|_2 \quad (2.1) $$

where $\|A\|_2$ denotes the spectral norm, i.e., largest singular value, of A. The Frobenius norm of A is denoted by $\|A\|_F$. If A is a $n \times n$ matrix with real eigenvalues, we order them as $\lambda_n(A) \le \dots \le \lambda_1(A)$. Positive (absolute) constants are denoted by c, C, $\tilde{C}$ etc. with suitable choice of indices. For $a, b \ge 0$, we write $a \lesssim b$ if there exists a constant $C > 0$ such that $a \le Cb$. Moreover, we write $a \asymp b$ if $a \lesssim b$ and $b \lesssim a$.

## 2.1 BTL model in a dynamic setting

Let us formally introduce our model for dynamic pairwise comparisons, inspired by the Bradley-Terry-Luce (BTL) model. We consider a set of items $[n] = \{1, 2, \dots, n\}$, with a certain quality at each time $t \in [0, 1]$, represented by the weight vector $w_t^* = (w_{t,1}^*, \dots, w_{t,n}^*)^{\top} \in \mathbb{R}^n$ with $w_{t,i}^* > 0$ for each $i \in [n]$. Let us denote

$$ b(t) := \max_{i,j \in [n]} \frac{w_{t,i}^*}{w_{t,j}^*} \quad \text{for all } t \in [0, 1] \qquad (2.2) $$

where we will require that $b(t)$ is finite for each $t \in [0, 1]$.

Our data consists of pairwise comparisons on this set of items at times $t \in \mathcal{T} = \{\frac{i}{T} | i = 0, \dots, T\}$. The outcomes at each time $t \in \mathcal{T}$ are gathered into an undirected comparison graph $G_t = ([n], \mathcal{E}_t)$ where $\mathcal{E}_t$ is the set of edges. While the set of items $[n]$ is supposed to be the same throughout, the compared items, i.e. the set of edges $\mathcal{E}_t$ can change with time.

To model such data, we use the BTL model at each time $t \in \mathcal{T}$. This model posits that the probability that an item $j$ wins over an item $i$ is proportional to its strength. At each time $t \in \mathcal{T}$, for each pair of compared items $\{i, j\} \in \mathcal{E}_t$, we perform $L$ independent comparisons. Their outcomes are independent Bernoulli variables, defined for $l \in \{1, \dots, L\}$ by $y_{ij}^{(l)}(t)$ where

$$ \mathbb{P}(y_{ij}^{(l)}(t) = 1) = \frac{w_{t,j}^{*}}{w_{t,i}^{*} + w_{t,j}^{*}}. $$

The proportion of times $j$ won over $i$ at time $t$ is given by $y_{ij}(t) = \frac{1}{L} \sum_{l=1}^{L} y_{ij}^{(l)}(t)$ and the corresponding true proportion is denoted by

$$ y_{ij}^{*}(t) := \mathbb{E}[y_{ij}(t)] = \frac{w_{t,j}^{*}}{w_{t,i}^{*} + w_{t,j}^{*}} \quad \forall i \neq j. $$

**Smooth evolution of pairwise outcomes.** Our goal is to recover $w_t^*$ at any time $t \in [0, 1]$. Suppose for the moment that $t$ is on the grid, then if $G_t$ is connected, $w_t^*$ is identifiable up to a positive scaling. However in our dynamic setting $G_t$ can be very sparse, and is not necessarily connected. Therefore for meaningful recovery of $w_t^*$, we need to make additional assumptions on the evolution of the pairwise outcomes over time. To this end, we make the following smoothness assumption.

**Assumption 1 (Lipschitz smoothness).** There exists $M \ge 0$ such that

$$ |y_{ij}^{*}(t) - y_{ij}^{*}(t')| \le M|t - t'| \quad \forall t, t' \in [0, 1], \quad i \neq j \in [n]. \qquad (2.3) $$
This assumption suggests that the pairwise outcomes at nearby time instants are similar, hence it is plausible that $w_t^*$ (at any $t \in [0, 1]$) could be estimated by utilizing the data lying in a neighborhood of $t$. To formalize this intuition, let us define a neighborhood at any time $t$ by

$$ \mathcal{N}_{\delta}(t) := \{t' \in \mathcal{T} \mid |t - t'| \le \frac{\delta}{T}\}. \qquad (2.4) $$

where $\delta \in [0, T]$. Note that if $\delta < \frac{1}{2}$, then there exists some values of $t$ for which $\mathcal{N}_{\delta}(t)$ is empty; hence we will consider $\delta \in [1/2, T]$. It is easy to verify that $\delta \le |\mathcal{N}_{\delta}(t)| \le 4\delta$ for all $t \in [0, 1]$. It will also be useful to denote

$$ b_{\max,\delta}(t) = \max_{t' \in \mathcal{N}_{\delta}(t)} b(t'). \qquad (2.5) $$

**Neighborhood graph.** For any time $t \in [0, 1]$, the data contained in the neighborhood $\mathcal{N}_{\delta}(t)$ can be gathered into a union graph, defined as $G_{\delta}(t) = ([n], \mathcal{E}_{\delta}(t))$ where $\mathcal{E}_{\delta}(t) = \cup_{t' \in \mathcal{N}_{\delta}(t)} \mathcal{E}_{t'}$. The maximum and minimum degree of a vertex in $G_{\delta}(t)$ will be denoted by $d_{\delta,\max}(t)$ and $d_{\delta,\min}(t)$ respectively. For each edge $\{i, j\} \in \mathcal{E}_{\delta}(t)$, it will also be useful to denote the time instants in $\mathcal{N}_{\delta}(t)$ where $i$ and $j$ are compared as

$$ \mathcal{N}_{ij,\delta}(t) = \{t' \in \mathcal{N}_{\delta}(t) \mid \{i, j\} \in \mathcal{E}_{t'}\}, $$

along with the quantities

$$ N_{\max,\delta}(t) = \max_{\{i,j\} \in \mathcal{E}_{\delta}(t)} |\mathcal{N}_{ij,\delta}(t)| \quad \text{and} \quad N_{\min,\delta}(t) = \min_{\{i,j\} \in \mathcal{E}_{\delta}(t)} |\mathcal{N}_{ij,\delta}(t)|. $$

**General recovery idea.** Given the above setup, a general idea for recovering $w_t^*$ is to first form the union graph $G_{\delta}(t)$, and to then compute the statistics

$$ \bar{y}_{ij}(t) := \frac{1}{|\mathcal{N}_{ij,\delta}(t)|} \sum_{t' \in \mathcal{N}_{ij,\delta}(t)} y_{ij}(t') \quad \text{for } \{i,j\} \in \mathcal{E}_{\delta}(t). \qquad (2.6) $$

In the limit $L \to \infty$, we have for each $\{i, j\} \in \mathcal{E}_{\delta}(t)$ that

$$ \bar{y}_{ij}(t) \xrightarrow{L \to \infty} \bar{y}_{ij}^*(t) = \frac{1}{|\mathcal{N}_{ij,\delta}(t)|} \sum_{t' \in \mathcal{N}_{ij,\delta}(t)} y_{ij}^*(t'). $$

Due to Assumption 1 we know that

$$ |\bar{y}_{ij}^*(t) - y_{ij}^*(t)| \le M\delta/T, \quad \text{for } \{i, j\} \in \mathcal{E}_{\delta}(t), $$

hence if $\delta = o(T)$, then for each $\{i, j\} \in \mathcal{E}_{\delta}(t)$, we have that $\bar{y}_{ij}^*(t)$ converges to $y_{ij}^*(t)$ as $T \to \infty$. Moreover, if the corresponding sequence of graphs $G_{\delta}(t)$ is connected, the identifiability of $w_t^*$ (up to a positive scaling) is ensured. The connectivity requirement on $G_{\delta}(t)$ is of course much weaker than requiring the individual graph(s) to be connected, which is also a key difference between the static and dynamic settings. Also note that the computed statistics $\bar{y}_{ij}(t)$ are nearest neighbor estimators as they average the data over a suitable neighborhood of $t$. More generally, one could also compute $\bar{y}_{ij}(t)$ via kernel smoothing as considered in [2].

With the above discussion in mind, a general scheme for recovering $w_t^*$ for any given $t \in [0, 1]$ would be to first form $G_{\delta}(t)$ for a suitable choice of the parameter $\delta$, then compute the statistics $\bar{y}_{ij}(t)$ as in (2.6), and finally apply any existing method for the static case using the comparison graph $G_{\delta}(t)$ and the data $(\bar{y}_{ij}(t))_{\{i,j\} \in \mathcal{E}_{\delta}(t)}$. We will focus on the Rank Centrality algorithm of Negahban et al. [23] – a popular spectral algorithm known to achieve state of the art performance – and adapt it to our dynamic setting.
## 2.2 Spectral dynamic ranking

The Rank Centrality [23] method is based on the connection between pairwise comparisons and a random walk on a directed graph. In the static case, a faithful estimation of the weight vector is given by the stationary distribution of the Markov chain induced by a suitably constructed transition matrix. In the dynamic setting, this method can be adapted for estimating $w_t^*$ by now constructing the transition matrix using $G_\delta(t)$ and the data $(\bar{y}_{ij}(t))_{\{i,j\} \in \mathcal{E}_\delta(t)}$, as in (2.6).

More precisely, we define a transition matrix $\hat{P}(t)$ on this graph using $(\bar{y}_{ij}(t))_{\{i,j\} \in \mathcal{E}_\delta(t)}$, with

$$ \hat{P}_{ij}(t) = \begin{cases} \frac{\bar{y}_{ij}(t)}{d_{\delta}(t)} = \frac{1}{d_{\delta}(t)} \left( \frac{1}{|\mathcal{N}_{ij,\delta}(t)|} \sum_{t' \in \mathcal{N}_{ij,\delta}(t)} y_{ij}(t') \right) & \text{if } \{i,j\} \in \mathcal{E}_{\delta}(t) \\ 1 - \frac{1}{d_{\delta}(t)} \left( \sum_{k \neq i} \frac{1}{|\mathcal{N}_{ik,\delta}(t)|} \sum_{t' \in \mathcal{N}_{ik,\delta}(t)} y_{ik}(t') \right) & \text{if } i = j \\ 0 & \text{otherwise} \end{cases} \quad (2.7) $$

where $d_\delta(t) \ge d_{\delta,\max}(t)$ is a suitably chosen normalization term. Then, one can easily verify that at each time $t$, $\hat{P}(t)$ is a transition matrix ($\hat{P}(t)$ is stochastic) corresponding to a Markov chain on a finite state space. Thus there always exists at least one stationary distribution associated to $\hat{P}(t)$. Moreover, stochastic matrices admit 1 as leading eigenvalue and so a candidate as the stationary distribution is its leading left eigenvector, i.e.

$$ \hat{\pi}(t)^{\top} = \hat{\pi}(t)^{\top} \hat{P}(t). $$

Besides, the vector of true weights $w_t^*$ we want to recover can be seen as the stationary distribution of a transition matrix on the union graph. Specifically, denoting $\pi^*(t) = \frac{w_t^*}{\sum_{i=1}^n w_{t,i}^*}$, one can easily show that $\pi^*(t)$ is the stationary distribution of the transition matrix

$$ \bar{P}_{ij}(t) = \begin{cases} \frac{1}{d_{\delta}(t)} \frac{w_{t,j}^*}{w_{t,i}^* + w_{t,j}^*} & \text{if } \{i,j\} \in \mathcal{E}_{\delta}(t) \\ 1 - \frac{1}{d_{\delta}(t)} \sum_{k \neq i} \frac{w_{t,k}^*}{w_{t,i}^* + w_{t,k}^*} & \text{if } i = j \\ 0 & \text{otherwise} \end{cases} \quad (2.8) $$

since $\bar{P}(t)$ and $\pi^*(t)$ verify the detailed balance equation of reversibility [13]

$$ \bar{P}_{ij}(t)\pi_i^*(t) = \bar{P}_{ji}(t)\pi_j^*(t) \quad \forall i, j \in [n]. $$

One can reasonably expect $\hat{\pi}(t)$ to be close to $\pi^*(t)$ as they are stationary distributions of $\hat{P}(t)$ and $\bar{P}(t)$ respectively, the latter of which are expected to be close. Indeed, one has the following bias-variance trade-off

$$ \hat{P}(t) - \bar{P}(t) = \underbrace{\hat{P}(t) - \mathbb{E}[\hat{P}(t)]}_{\text{variance}} + \underbrace{\mathbb{E}[\hat{P}(t)] - \bar{P}(t)}_{\text{bias}}. $$

where the variance term is typically expected to decrease with $\delta$ (due to averaging over $\mathcal{N}_\delta(t)$) while the bias term will scale as $O(\delta/T)$ (due to the smoothness assumption 1). Hence for a suitably chosen $\delta = o(T)$ we will then have (for $n, T$ large enough) $\hat{P}(t) \approx \bar{P}(t)$, which implies $\hat{\pi}(t) \approx \pi^*(t)$.
**Remark 1.** For meaningful recovery the vector $\pi^*(t)$ clearly has to be unique. This is the case if the associated Markov chain is irreducible which in turn is ensured by the connectivity of the underlying graph (here, the union graph $G_\delta(t)$) and the strict positivity of the weights on its edges [13]. The condition on the weights is guaranteed in our setup since $\tilde{P}_{ij}(t) > 0$ for each $\{i, j\} \in \mathcal{E}_\delta(t)$ (indeed, $w_{t,i}^* > 0$ for each $t \in [0, 1]$ and $i \in [n]$).

Based on the above discussion we can outline the steps of our method for ranking in the dynamic setting in the form of Algorithm 1. Our goal now is to establish conditions under which $\hat{\pi}(t)$ is

<caption>**Algorithm 1** Spectral algorithm for dynamic ranking (Dynamic Rank Centrality)</caption>

1: **Input:** Grid $\mathcal{T} \subset [0, 1]$, and a given time $t \in [0, 1]$. For each $t' \in \mathcal{T}$: comparison graph $G_{t'}$, results of comparisons as statistics $(y_{ij}(t'))_{\{i,j\} \in \mathcal{E}_{t'}}$.

2: Form the neighbourhood graph $G_\delta(t) = ([n], \mathcal{E}_\delta(t))$ where $\mathcal{E}_\delta(t) = \cup_{t' \in \mathcal{N}_\delta(t)} \mathcal{E}_{t'}$, and $\mathcal{N}_\delta(t)$ is as in (2.4).

3: Compute the transition matrix $\hat{P}(t)$ as in (2.7) with $d_\delta(t) \ge d_{\delta,\max}(t)$.

4: Compute the leading left eigenvector $\hat{\pi}(t)$ of $\hat{P}(t)$.

5: **Output:** $\hat{\pi}(t) \in \mathbb{R}^n$.

close to $\pi^*(t)$ under the $l_2$ and $l_\infty$ norms. These results are summarized in the next section. In particular, we will strive to establish consistency results (i.e., the error approaching zero) when the grid size $T \to \infty$.

Before proceeding, we need to define some additional quantities related to the union graph $G_\delta(t)$ which will appear in the following sections. Let $\mathcal{L}_\delta(t) = D_\delta^{-1}(t)A_\delta(t)$ denote the random walk Laplacian of $G_\delta(t)$, where $D_\delta(t)$ is the diagonal matrix of vertex degrees, and $A_\delta(t)$ is its adjacency matrix. We denote $\xi_\delta(t) = 1 - \lambda_{\max}(\mathcal{L}_\delta(t))$ where

$$ \lambda_{\max}(\mathcal{L}_\delta(t)) := \max \{\lambda_2(\mathcal{L}_\delta(t)), -\lambda_n(\mathcal{L}_\delta(t))\}, $$

is the second largest eigenvalue (in absolute value) of $\mathcal{L}_\delta(t)$. Note that $\mathcal{L}_\delta(t)$ has real eigenvalues since it is similar to the symmetric Laplacian $D_\delta^{-1/2}(t)A_\delta(t)D_\delta^{-1/2}(t)$.

# 3 Main results

In Section 3.1, we present bounds on the $l_2$ error $\|\hat{\pi}(t) - \pi^*(t)\|_2$, while Section 3.2 contains our bounds on the $l_\infty$ error $\|\hat{\pi}(t) - \pi^*(t)\|_\infty$.

## 3.1 $l_2$ error bound

For a given sequence of graphs $(G(t'))_{t' \in \mathcal{T}}$, the following theorem provides an explicit $l_2$ error bound (holding w.h.p) which in particular highlights the dependence on parameters related to the union graph $G_\delta(t)$, the grid size $T$ and the neighborhood size $\delta$.

**Theorem 1.** For any given $t \in [0, 1]$, suppose that $\delta \in [\frac{1}{2}, T]$ is such that $n \ge c_1 \log n$ and $\xi_\delta(t) > 0$ for some constant $c_1 > 0$. Choosing $d_\delta(t) \ge d_{\delta,\max}(t)$, there exist constants $\tilde{C}_1 \ge 15, \tilde{C}_2 \ge 1$ such that if

$$ \tilde{C}_1 \sqrt{\frac{N_{\max,\delta}(t) d_{\delta,\max}(t) \log n}{L d_\delta^2(t) N_{\min,\delta}^2(t)}} + 4 \frac{M \delta |\mathcal{E}_\delta(t)|}{T d_{\delta,\max}(t)} \le \frac{\xi_\delta(t) d_{\delta,\min}(t)}{8 d_\delta(t) b^{7/2}(t)}, \quad (3.1) $$
then it holds with probability at least $1 - O(n^{-10})$ that

$$ \frac{\|\hat{\pi}(t) - \pi^*(t)\|_2}{\|\pi^*(t)\|_2} \le 32 \frac{M\delta|\mathcal{E}_\delta(t)|b^{7/2}(t)d_\delta(t)}{T\xi_\delta(t)d_{\delta,\min}(t)d_{\delta,\max}(t)} + 8\tilde{C}_2 \frac{b^{9/2}(t)}{\xi_\delta(t)d_{\delta,\min}(t)} \sqrt{\frac{N_{\max,\delta}(t)d_{\delta,\max}(t)}{LN_{\min,\delta}^2(t)}}. \quad (3.2) $$

Let us make the following observations.

1. The first term in the RHS of (3.2) corresponds to the bias and arises from the regularity assumption 1, while the second term therein is the variance term. Moreover, note that the error depends on $\delta$ – either explicitly, or through certain quantities such as $d_{\delta,\max}(t)$, $N_{\min,\delta}(t)$ etc. In order to obtain a more explicit dependence in terms of $\delta$, we will need to make specific assumptions on the graphs $G_{t'}$, $t' \in \mathcal{T}$. Below, we will consider the setting where the graphs are Erdös-Renyi graphs and derive explicit conditions on $\delta$ that lead to consistency with respect to $\mathcal{T}$.

2. In the static case, we have $t = t'$ for some $t' \in \mathcal{T}$ and only the graph $G_{t'}$ is observed. Then $M = 0$ and $\delta = 1/2$, so $N_{\min,\delta}(t)$, $N_{\max,\delta}(t) \equiv 1$. Denoting $d_{\min}(t)$, $d_{\max}(t)$, $\xi(t)$ to be the corresponding quantities with the $\delta$ suffix suppressed, condition (3.1) is satisfied for $L$ large enough. Moreover, the error bound is then $O(\frac{b^{9/2}(t)}{\xi(t)d_{\min}(t)\sqrt{\frac{d_{\max}(t)}{L}}})$ which matches the $\ell_2$ bound of Negahban et al. [23, Theorem 1] with the $\sqrt{\log n}$ factor therein removed, but with an extra $b^2(t)$ factor.

Now we consider the important case where the comparison graphs are Erdös-Renyi graphs, i.e., $G_{t'} = \mathcal{G}(n, p(t'))$ for each $t' \in \mathcal{T}$. It is easily seen that the union graph $G_\delta(t)$ is also then Erdös-Renyi denoted by $\mathcal{G}(n, p_\delta(t))$ where $p_\delta(t)$ is given by

$$ p_{\delta}(t) = 1 - \prod_{t' \in \mathcal{N}_{\delta}(t)} (1 - p(t')). \quad (3.3) $$

In this setting, the bound in Theorem 1 can be simplified using concentration results for parameters related to $G_\delta(t)$ (see Lemma 11). Specifically, we have that if $p_\delta(t) \gtrsim \log n/n$, then w.h.p

$$ \frac{np_\delta(t)}{2} \le d_{\delta,\max}(t), d_{\delta,\min}(t) \le \frac{3np_\delta(t)}{2}; \quad \xi_\delta(t) \ge \frac{1}{2}; \quad |\mathcal{E}_\delta(t)| \le 2n^2p_\delta(t). $$

In particular, we will choose the normalization factor $d_\delta(t) = 3np_\delta(t)$ which is a valid choice (w.h.p). Lemma 11 also states that if $p_{\delta,\text{sum}}(t) \gtrsim \log n$, then $N_{\max,\delta}(t)$, $N_{\min,\delta}(t) \asymp p_{\delta,\text{sum}}(t)$ w.h.p. These considerations lead to the following simplification of Theorem 1.

**Theorem 2.** Suppose that $G_{t'} \sim \mathcal{G}(n, p(t'))$ for all $t' \in \mathcal{T}$ so that $G_\delta(t) \sim \mathcal{G}(n, p_\delta(t))$ (for any given $t \in [0, 1]$) with $p_\delta(t)$ as in (3.3), and denote $p_{\delta,\text{sum}}(t) := \sum_{t' \in \mathcal{N}_\delta(t)} p(t')$. Choosing $d_\delta(t) = 3np_\delta(t)$, let $\delta \in [\frac{1}{2}, T]$ be such that $n \ge c_1 \log n$, $np_\delta(t) \ge c_0 \log n$, and $p_{\delta,\text{sum}}(t) \ge c_2 \log n$ with constant $c_1 > 0$ as in Theorem 1, and constants $c_0, c_2 \ge 1$. Then for constants $\tilde{C}_1, \tilde{C}_2$ as in Theorem 1, if

$$ 2\tilde{C}_1 \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} + 16 \frac{M\delta n}{T} \le \frac{1}{96b^{7/2}(t)} \quad (3.4) $$

holds, we have with probability at least $1 - O(n^{-10})$ that

$$ \frac{\|\hat{\pi}(t) - \pi^*(t)\|_2}{\|\pi^*(t)\|_2} \le 1536 \frac{M\delta nb^{7/2}(t)}{T} + 64\tilde{C}_2 b^{9/2}(t) \sqrt{\frac{3}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}}. $$
The following remarks are in order.

1. As can be seen, the bias term is $O(\frac{n\delta}{T})$ while the variance term scales as $O(\frac{1}{\sqrt{Lnp_{\delta}(t)p_{\delta,\text{sum}}(t)}})$. Hence if $p_{\delta,\text{sum}}(t)$ grows with $\delta$, then the variance error will reduce as $\delta$ increases.

2. In the static case we observe a single comparison graph $G_{t'}$ (for $t' \in T$) with $t = t'$. Thus $p_{\delta,\text{sum}}(t) \equiv 1$ and the condition $p_{\delta,\text{sum}}(t) \gtrsim \log n$ is not needed, while $p_{\delta}(t) = p(t)$. Hence, if $p(t) \gtrsim \frac{\log n}{n}$ and $L$ is suitably large, the $\ell_2$ error is bounded by $O(\frac{1}{\sqrt{Lnp(t)}})$, which corresponds to the bound obtained by Chen et al. [7, Theorem 9]. So our result is coherent with existing results for the static case for Erdös-Renyi graphs.

3. The choice $d_{\delta}(t) = 3np_{\delta}(t)$ ensures that $d_{\delta,\text{max}}(t) \le d_{\delta}(t)$ w.h.p. In fact, we could have chosen $d_{\delta}(t)$ to be a constant ($\ge 1$) multiple of $d_{\delta,\text{max}}(t)$ as well. However for the $\ell_{\infty}$ analysis later on, it will be crucial to choose $d_{\delta}(t)$ as a constant times $np_{\delta}(t)$ for technical reasons arising in the analysis. Similar considerations for the choice of the normalization factor exist in the static setting as well (see [7, 6]). Note that this choice of $d_{\delta}(t)$ requires us to know $p_{\delta}(t)$, but in case we don't know $p_{\delta}(t)$ in practice, we can instead use its empirical estimate which can be easily computed.

We now derive an appropriate choice for $\delta$ that leads to an $\ell_2$ error rate of $O(T^{-1/3})$. To this end, we first need to explicitly show the dependence on $\delta$ for $p_{\delta}(t), p_{\delta,\text{sum}}(t)$. Let us assume for simplicity that

$$p_{\min} := \min_{t' \in T} p(t') > 0. \qquad (3.5)$$

Since $\delta \le |\mathcal{N}_{\delta}(t)| \le 4\delta$, we have for all $t \in [0, 1]$ that $p_{\delta,\text{sum}}(t) \ge \delta p_{\min}$. Besides, as shown in Proposition 1, $p_{\delta}(t) \gtrsim \min\{1, \delta p_{\min}\}$. Hence if $\delta p_{\min} \gtrsim \log n$ then it implies

$$p_{\delta,\text{sum}}(t) \gtrsim \log n \quad \text{and} \quad p_{\delta}(t) \gtrsim 1 (\ge \log n/n)$$

meaning that the conditions on $p_{\delta}(t)$ and $p_{\delta,\text{sum}}(t)$ in Theorem 2 are satisfied.

**Remark 2.** The condition $p_{\delta,\text{sum}}(t) \gtrsim \log n$ is only needed to ensure that $N_{\min,\delta}(t), N_{\max,\delta}(t)$ concentrate around $p_{\delta,\text{sum}}(t)$. In fact, if we replace $|\mathcal{N}_{ij,\delta}(t)|$ with $p_{\delta,\text{sum}}(t)$ in (2.7), then this condition would disappear and we would only be left with the requirement $p_{\delta}(t) \gtrsim \log n/n$. This would of course impose a less stringent lower bound on $T$ as it would then suffice that $\delta p_{\min} \gtrsim \frac{\log n}{n}$. Since our setup is not limited to Erdös-Renyi graphs, we prefer to keep $|\mathcal{N}_{ij,\delta}(t)|$ in (2.7).

**Corollary 1.** Under the same notations as in Theorem 2, for all $t \in [0, 1]$ suppose that $n \gtrsim \log n$ and $p_{\min}$ is as in (3.5). Choosing $\delta = \min\left\{\frac{(b(t))^{\frac{2}{3}}}{(2M)^{\frac{2}{3}}n(Lp_{\min})^{\frac{1}{3}}}T^{2/3}, T\right\}$ and $d_{\delta}(t) = 3np_{\delta}(t)$, if $T$ is such that $\delta \gtrsim \frac{\log n}{p_{\min}}$ and

$$\left(\sqrt{\frac{\log n}{Lnp_{\min}}}\right) \max \left\{ \frac{(2M)^{1/3}n^{1/2}(Lp_{\min})^{1/6}}{b^{1/3}(t)T^{1/3}}, \frac{1}{\sqrt{T}} \right\} + Mn \min \left\{ \frac{b^{2/3}(t)}{(2M)^{2/3}n(Lp_{\min})^{1/3}T^{1/3}}, 1 \right\} \lesssim \frac{1}{b^{7/2}(t)}, \\ \text{then with probability at least } 1 - O(n^{-10}),$$

$$\frac{\|\hat{\pi}(t) - \pi^*(t)\|_2}{\|\pi^*(t)\|_2} \lesssim Mnb^{7/2}(t) \min \left\{ \frac{b^{2/3}(t)}{(2M)^{2/3}n(Lp_{\min})^{1/3}T^{1/3}}, 1 \right\} \\ + \frac{b^{9/2}(t)}{\sqrt{Lnp_{\min}}} \max \left\{ \frac{(2M)^{1/3}n^{1/2}(Lp_{\min})^{1/6}}{b^{1/3}(t)T^{1/3}}, \frac{1}{\sqrt{T}} \right\}.$$
The following observations are useful to note.

1. When $M > 0$, Corollary 1 states that for $\delta = \Theta(T^{2/3})$, if $n, T$ are large enough, then w.h.p $\|\hat{\pi}(t) - \pi^*(t)\|_2 = O(T^{-1/3})$. This matches the rate for the pointwise risk for estimating univariate Lipschitz functions (see for e.g. [24, Theorem 1.3.1]).

2. If $M = 0$ then $\delta = T$ which makes sense since $y_{ij}^*(t)$ is a constant function for each $i \neq j$. Indeed, the problem is then the same as the setting where the comparison graph is $\cup_{t' \in \mathcal{T}} G_{t'}$, and we observe (a potentially different number of) i.i.d pairwise outcomes for each given edge in this graph. In this case, the corollary states that provided $T$ is large enough, the $\ell_2$ error is $\lesssim \frac{b^{9/2}(t)}{\sqrt{L np_{\min} T}}$. This is logically faster than the $T^{-1/3}$ nonparametric rate, and is analogous to the optimal $\ell_2$ bound for Erdös Renyi graphs in the static setting (see [7, Theorem 9]).

## 3.2 $\ell_\infty$ error bound

We now discuss our results for bounding the $\ell_\infty$ error $\|\hat{\pi}(t) - \pi^*(t)\|_\infty$ at any given time $t$. Such bounds are particularly desirable in the context of ranking as they lead to guarantees for recovering the ranks of the items. We will assume that all the comparison graphs (at each $t' \in \mathcal{T}$) are Erdös-Renyi graphs. The following theorem is the $\ell_\infty$ counterpart of Theorem 2, its proof is outlined in Section 5.

**Theorem 3.** Under the notation and assumptions of Theorem 2, there exists a constant $\tilde{C}_3 \ge 1$ such that if additionally

$$96b^{\frac{5}{2}}(t) \left( \frac{4M\delta}{T} + \tilde{C}_3 \sqrt{\frac{\log n}{np_\delta(t)}} \right) \le \frac{1}{2}, \quad (3.6)$$

then there exist constants $\tilde{C}_4, \tilde{C}_5, \tilde{C}_6 \ge 1$ such that with probability at least $1 - O(n^{-9})$,

$$\frac{\|\hat{\pi}(t) - \pi^*(t)\|_{\infty}}{\|\pi^*(t)\|_{\infty}} \le \frac{12b_{\max,\delta}(t)}{1 - \tilde{C}_4 b_{\max,\delta}(t) \sqrt{\frac{\log n}{np_\delta(t)}}} \left( \tilde{C}_5 \gamma_{n,\delta}(t) \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} + \tilde{C}_6 \frac{Mn\delta b^{\frac{7}{2}}(t)}{T} \right),$$

where $\gamma_{n,\delta}(t) := (1 + \frac{b^{\frac{5}{2}}(t)}{\sqrt{\log n}} \max\{b^2(t), \frac{\log n}{\sqrt{np_\delta(t)}}\})$ and $b_{\max,\delta}(t)$ is as in (2.5).

As before for Theorem 2, let us interpret Theorem 3 for the static setting where $t = t'$ for some $t' \in \mathcal{T}$, and only $G_{t'}$ is observed. Then Theorem 3 states that if $np(t) \gtrsim b^5(t) \log n$, and $n, L$ are large enough, then w.h.p, the $\ell_\infty$ error is

$$\frac{\|\hat{\pi}(t) - \pi^*(t)\|_{\infty}}{\|\pi^*(t)\|_{\infty}} \lesssim b(t) \left( 1 + \frac{b^{\frac{5}{2}}(t)}{\sqrt{\log n}} \max\{b^2(t), \frac{\log n}{\sqrt{np_\delta(t)}}\} \right) \sqrt{\frac{\log n}{Lnp(t)}}.$$

Hence if $b(t) = O(1)$ then the bound is $O(\sqrt{\frac{\log n}{Lnp(t)}})$ which matches the corresponding bound of Chen et al. [7, Theorem 3].

Let us now denote

$$b_{\max} := \max_{t' \in \mathcal{T}} b(t'), \quad \gamma_n(t) := \left( 1 + \frac{b^{\frac{5}{2}}(t)}{\sqrt{\log n}} \max\{b^2(t), \frac{\log n}{\sqrt{n}}\} \right) \quad (3.7)$$

so that $b_{\max,\delta}(t) \le b_{\max}$ and $\gamma_{n,\delta}(t) \lesssim \gamma_n(t)$. Then as for the $\ell_2$ case, we can derive a value for $\delta$ that leads to a $\ell_\infty$ error rate of $T^{-\frac{1}{3}}$.
**Corollary 2.** Under the same notations as in Theorem 3, for all $t \in [0, 1]$ suppose that $n \gtrsim \log n$, $p_{\min}$ is as in (3.5) and $b_{\max}$, $\gamma_n(t)$ are as in (3.7). Choosing $\delta = \min\left\{\frac{(\gamma_n(t))^{\frac{2}{3}}(\log n)^{\frac{1}{3}}}{(2M)^{\frac{2}{3}}nb^{\frac{7}{3}}(t)(Lp_{\min})^{\frac{1}{3}}}T^{2/3}, T\right\}$ and $d_{\delta}(t) = 3np_{\delta}(t)$, if $T$ is such that $\delta \gtrsim \frac{\log n}{p_{\min}}$ and

$$ \min \left\{ \left( \frac{(\gamma_n(t))^{\frac{2}{3}}(M \log n)^{\frac{1}{3}}}{2^{\frac{2}{3}}nb^{\frac{7}{3}}(t)(Lp_{\min})^{\frac{1}{3}}} \right) \frac{1}{T^{1/3}}, M \right\} + \sqrt{\frac{\log n}{n}} \lesssim \frac{1}{b^{5/2}(t)}, $$

then with probability at least $1 - O(n^{-9})$,

$$ \frac{\|\hat{\pi}(t) - \pi^*(t)\|_{\infty}}{\|\pi^*(t)\|_{\infty}} \lesssim \left( \frac{b_{\max}}{1 - b_{\max}\sqrt{\frac{\log n}{n}}} \right) \left[ \gamma_n(t) \sqrt{\frac{\log n}{Lnp_{\min}}} \max \left\{ \frac{(2M)^{1/3}n^{1/2}b^{7/6}(t)(Lp_{\min})^{1/6}}{\gamma_n^{1/3}(t)(\log n)^{1/6}T^{1/3}}, \frac{1}{\sqrt{T}} \right\} \right. \\ \left. + Mnb^{7/2}(t) \min \left\{ \frac{\gamma_n^{2/3}(t)(\log n)^{1/3}}{(2M)^{2/3}nb^{7/3}(t)(Lp_{\min})^{1/3}T^{1/3}}, 1 \right\} \right]. $$

Let us make the following useful observations.

1. If $b(t') = O(1)$ for some $t' \in \mathcal{T}$ then it is not difficult to verify that this implies $b(t') = O(1)$ for each $t' \in \mathcal{T}$ (and hence $b_{\max} = O(1)$), due to Assumption 1. Additionally, this implies that $\gamma_n(t) = O(1)$.

2. When $M > 0$ and $b_{\max} = O(1)$, Corollary 2 asserts that for $\delta = \Theta(T^{2/3})$, if $n, T$ are large enough, then w.h.p $\|\hat{\pi}(t) - \pi^*(t)\|_{\infty} = O(T^{-1/3})$. This matches the rate for the pointwise risk for estimating univariate Lipschitz functions (see for e.g. [24, Theorem 1.3.1]).

# 4 $\ell_2$-analysis of the spectral estimator

We now describe the main ideas that lead to the $\ell_2$ bound in Theorem 1. We will essentially proceed in three steps following the ideas in [7].

$$ \begin{aligned} \|\hat{\pi}(t) - \pi^*(t)\|_2 &\stackrel{(i)}{\le} \frac{1}{\sqrt{\pi_{\min}^*(t)}} \|\hat{\pi}(t) - \pi^*(t)\|_{\pi^*(t)} \\ &\stackrel{(ii)}{\le} \frac{8d_{\delta}(t)b^{7/2}(t)}{\xi_{\delta}(t)d_{\delta,\min}(t)} \|\pi^*(t)^{\top}(\hat{P}(t) - \bar{P}(t))\|_2 \end{aligned} \quad (4.1) $$

$$ \stackrel{(iii)}{\le} \frac{8d_{\delta}(t)b^{7/2}(t)}{\xi_{\delta}(t)d_{\delta,\min}(t)} \left( 4 \frac{M\delta|\mathcal{E}_{\delta}(t)|}{Td_{\delta,\max}(t)} + \tilde{C}_2 \sqrt{\frac{N_{\max,\delta}(t)d_{\delta,\max}(t)b^2(t)}{Ld_{\delta}^2(t)N_{\min,\delta}^2(t)}} \right) \|\pi^*(t)\|_2. \quad (4.2) $$

(i) The first step is easy to verify, due to the definition of the norm $\|\cdot\|_{\pi^*(t)}$

$$ \|\pi^*(t) - \hat{\pi}(t)\|_{\pi^*(t)}^2 = \sum_{i=1}^{n} \pi_i^*(t) (\pi_i^*(t) - \hat{\pi}_i(t))^2 \ge \pi_{\min}^*(t) \|\pi^*(t) - \hat{\pi}(t)\|_2^2. $$

(ii) This is given by the combination of Lemmas 1, 2 and 3 which in turn are derived using [7, Theorem 8] and [23, Lemma 6].

(iii) For this step, we can decompose $\hat{P}(t) - \bar{P}(t) = \Delta(t) + \Delta_1(t)$ as in (4.5). Then one has to bound $\|\pi^*(t)^{\top}\Delta(t)\|_2$ and $\|\pi^*(t)^{\top}\Delta_1(t)\|_2$. The second term is completely deterministic and can be bounded using Assumption 1. A bound on $\|\pi^*(t)^{\top}\Delta(t)\|_2$ is found following the same steps as in the proof of [6, Theorem 9].
## 4.1 Proof of Theorem 1

A bound on $||\pi^*(t) - \hat{\pi}(t)||_{\pi^*(t)}$ is given by [7, Theorem 8], which is recalled in Appendix A. Denoting $\lambda_{\max}(\bar{P}(t))$ to be the second largest eigenvalue of $\bar{P}(t)$ in absolute value, i.e.,

$$ \lambda_{\max}(\bar{P}(t)) = \max \{ \lambda_2(\bar{P}(t)), -\lambda_n(\bar{P}(t)) \}, $$

this theorem gives the bound

$$ ||\hat{\pi}(t) - \pi^*(t)||_{\pi^*(t)} \le \frac{||\pi^*(t)^{\top}(\hat{P}(t) - \bar{P}(t))||_2}{1 - \lambda_{\max}(\bar{P}(t)) - ||\hat{P}(t) - \bar{P}(t)||_{\pi^*(t)}} \quad (4.3) $$

provided that the following condition holds.

$$ ||\hat{P}(t) - \bar{P}(t)||_{\pi^*(t)} < 1 - \lambda_{\max}(\bar{P}(t)). \quad (4.4) $$

First let us note that these eigenvalues are real, and so (4.4) is well defined. Indeed, denoting $\Pi^*(t) = \text{diag}(\pi^*(t))$ and $S = \Pi^*(t)^{1/2}\bar{P}(t)\Pi^*(t)^{-1/2}$, $S$ is similar to $\bar{P}(t)$, and $S$ is symmetric due to the reversibility of $\bar{P}(t)$.

To prove (4.4), we will use results similar to [23, Lemma's 3,4]. The main idea is to write the following decomposition

$$ \hat{P}(t) = \bar{P}(t) + \underbrace{\hat{P}(t) - \hat{P}^*(t)}_{\Delta(t)} + \underbrace{\hat{P}^*(t) - \bar{P}(t)}_{\Delta_1(t)} \quad (4.5) $$

where $\hat{P}^*(t) = \mathbb{E}\hat{P}(t)$ whose entries are given by

$$ \hat{P}_{ij}^*(t) = \begin{cases} \frac{1}{d_{\delta}(t)} \left( \frac{1}{|\mathcal{N}_{ij,\delta}(t)|} \sum_{t' \in \mathcal{N}_{ij,\delta}(t)} \frac{w_{t',j}^*}{w_{t',i}^* + w_{t',j}^*} \right) & \text{if } \{i,j\} \in \mathcal{E}_{\delta}(t) \\ 1 - \frac{1}{d_{\delta}(t)} \sum_{k \neq i} \frac{1}{|\mathcal{N}_{ik,\delta}(t)|} \sum_{t' \in \mathcal{N}_{ik,\delta}(t)} \frac{w_{t',k}^*}{w_{t',i}^* + w_{t',k}^*} & \text{if } i = j \\ 0 & \text{otherwise.} \end{cases} $$

We now provide bounds on $||\Delta_1(t)||_2$ and $||\Delta(t)||_2$ in Lemma's 1 and 2 respectively. The proofs of all results in this section are outlined in Appendix C.

**Lemma 1** (Bound on $||\Delta_1(t)||_2$). It holds that $||\Delta_1(t)||_2 \le 4 \frac{M\delta|\mathcal{E}_{\delta}(t)|}{Td_{\delta,\max}(t)}$.

**Lemma 2** (Bound on $||\Delta(t)||_2$). There exists a constant $\tilde{C}_1 \ge 15$ such that

$$ ||\Delta(t)||_2 \le \tilde{C}_1 \sqrt{\frac{N_{\max,\delta}(t)d_{\delta,\max}(t) \log n}{Ld_{\delta}^2(t)N_{\min,\delta}^2(t)}} $$

with probability at least $1 - O(n^{-10})$.

The proof of Lemma 1 follows from the smoothness condition in Assumption 1, while the proof of Lemma 2 follows the proof steps of [23, Lemma 3]. Next, we show that if $\xi_{\delta}(t) > 0$ (which implies that $G_{\delta}(t)$ is connected) and if the perturbation $||\hat{P}(t) - \bar{P}(t)||_2$ is sufficiently small, then we can ensure (4.4).
**Lemma 3.** Denoting $\rho(t) = \lambda_{\max}(\bar{P}(t)) + \|\hat{P}(t) - \bar{P}(t)\|_2\sqrt{\frac{\pi_{\max}^*(t)}{\pi_{\min}^*(t)}}$, recall that $\xi_{\delta}(t) = \lambda_{\max}(\mathcal{L}_{\delta}(t))$ (where $\mathcal{L}_{\delta}(t)$ is the Laplacian of $G_{\delta}(t)$) and $b(t) := \max_{i,j} \frac{w_{t,i}^*}{w_{t,j}^*} = \frac{\pi_{\max}^*(t)}{\pi_{\min}^*(t)}$. If $\xi_{\delta}(t) > 0$ and $d_{\delta}(t) \ge d_{\delta,\max}(t)$, then we have that $1 - \lambda_{\max}(\bar{P}(t)) \ge \frac{\xi_{\delta}(t)d_{\delta,\min}(t)}{4d_{\delta}(t)b^3(t)}$. Moreover, if

$$ \|\hat{P}(t) - \bar{P}(t)\|_2 \le \frac{\xi_{\delta}(t)d_{\delta,\min}(t)}{8d_{\delta}(t)b^{7/2}(t)} \quad (4.6) $$

then it holds that

$$ 1 - \rho(t) \ge \frac{\xi_{\delta}(t)d_{\delta,\min}(t)}{8d_{\delta}(t)b^3(t)} > 0. $$

The statement is analogous to that of [23, Lemma 4]. The bound on $1 - \lambda_{\max}(\bar{P}(t))$ is clearly the crucial statement, and requires using [23, Lemma 6]. We remark in passing that the dependence on $b(t)$ is $b^3(t)$ in Lemma 3, we could not verify the dependence stated in [23, Lemma 6] (which is $b^2(t)$). For completeness, we outline the proof of Lemma 3 in Appendix C.

Condition (4.6) is ensured via Lemma’s 1 and 2 (with high probability) whenever (3.1) holds. Then, (2.1) readily implies that

$$ 1 - \lambda_{\max}(\bar{P}(t)) - \|\hat{P}(t) - \bar{P}(t)\|_{\pi^*(t)} \ge 1 - \rho(t) \ge \frac{\xi_{\delta}(t)d_{\delta,\min}(t)}{8d_{\delta}(t)b^3(t)} > 0, \quad (4.7) $$

thus ensuring (4.4). Using (4.3) and (4.7) we finally obtain (4.1) as follows.

$$ \|\hat{\pi}(t) - \pi^*(t)\|_{\pi^*(t)} \le \frac{8d_{\delta}(t)b^3(t)}{\xi_{\delta}(t)d_{\delta,\min}(t)} \|\pi^*(t)^{\top}(\hat{P}(t) - \bar{P}(t))\|_{\pi^*(t)} \le \frac{8d_{\delta}(t)b^{7/2}(t)}{\xi_{\delta}(t)d_{\delta,\min}(t)} \|\pi^*(t)^{\top}(\hat{P}(t) - \bar{P}(t))\|_2 \quad (4.8) $$

where the last inequality uses (2.1).

Finally, we can bound $\|\pi^*(t)^{\top}(\hat{P}(t) - \bar{P}(t))\|_2$ using the decomposition in (4.5) along with the triangular inequality, leading to

$$ \|\pi^*(t)^{\top}(\hat{P}(t) - \bar{P}(t))\|_2 \le \|\pi^*(t)^{\top}\Delta(t)\|_2 + \|\pi^*(t)^{\top}\Delta_1(t)\|_2. $$

Bounds on $\|\pi^*(t)^{\top}\Delta(t)\|_2$ and $\|\pi^*(t)^{\top}\Delta_1(t)\|_2$ are provided in the following lemma's.

**Lemma 4.**

$$ \|\pi^*(t)^{\top}\Delta_1(t)\|_2 \le \|\pi^*(t)\|_2 \|\Delta_1(t)\|_2 \le 4 \frac{M\delta|\mathcal{E}_{\delta}(t)|}{Td_{\delta,\max}(t)} \|\pi^*(t)\|_2. \quad (4.9) $$

The statement follows directly from Lemma 1.

**Lemma 5.** There exist constants $c_1 > 0, \tilde{C}_2 \ge 1$ such that if $n \ge c_1 \log n$ then with probability at least $1 - O(n^{-10})$, we have that

$$ \|\pi^*(t)^{\top}\Delta(t)\|_2 \le \tilde{C}_2 \sqrt{\frac{N_{\max,\delta}(t)d_{\delta,\max}(t)b^2(t)}{Ld_{\delta}^2(t)N_{\min,\delta}^2(t)}} \|\pi^*(t)\|_2. $$

The proof of Lemma 5 follows the ideas in the proof of [7, Theorem 9]. Applying these bounds in (4.8) finally leads to the stated bound in Theorem 1.
## 4.2 Proof of Theorem 2

This theorem follows directly from Theorem 1 and from the properties of Erdös-Renyi graphs gathered in Lemma 11. Using (3.4) and Lemma 11 along with the choice $d_\delta(t) = 3np_\delta(t)$, it holds with probability at least $1 - O(n^{-10})$ that

$$ \tilde{C}_1 \sqrt{\frac{N_{\max,\delta}(t)d_{\delta,\max}(t) \log n}{Ld_\delta^2(t)N_{\min,\delta}^2(t)}} + 4 \frac{M\delta|\mathcal{E}_\delta(t)|}{Td_{\delta,\max}(t)} \le 2\tilde{C}_1 \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} + 16 \frac{M\delta n}{T} \le \frac{1}{96b^{\frac{7}{2}}(t)} \le \frac{\xi_\delta(t)d_{\delta,\min}(t)}{8d_\delta(t)b^{\frac{7}{2}}(t)}. $$

Hence, condition (3.1) is satisfied with high probability and Theorem 1 implies that

$$ \frac{\|\hat{\pi}(t) - \pi^*(t)\|_2}{\|\pi^*(t)\|_2} \le 32 \frac{M\delta|\mathcal{E}_\delta(t)|d_\delta(t)b^{7/2}(t)}{T\xi_\delta(t)d_{\delta,\min}(t)d_{\delta,\max}(t)} + 8\tilde{C}_2 \frac{b^{9/2}(t)}{\xi_\delta(t)d_{\delta,\min}(t)} \sqrt{\frac{N_{\max,\delta}(t)d_{\delta,\max}(t)}{LN_{\min,\delta}^2(t)}}. $$

Again, using Lemma 11, we can simplify the above bound so that with probability at least $1 - O(n^{-10})$,

$$ \frac{\|\hat{\pi}(t) - \pi^*(t)\|_2}{\|\pi^*(t)\|_2} \le 1536 \frac{M\delta n^2 p_\delta(t) b^{7/2}(t)}{T np_\delta(t)} + 64 \tilde{C}_2 \frac{b^{9/2}(t)}{np_\delta(t)} \sqrt{\frac{3p_{\delta,\text{sum}}(t)np_\delta(t)}{Lp_{\delta,\text{sum}}^2(t)}} \\ \le 1536 \frac{M\delta nb^{7/2}(t)}{T} + 64 \tilde{C}_2 b^{9/2}(t) \sqrt{\frac{3}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}}. $$

## 4.3 Proof of Corollary 1

Since $p_{\delta,\text{sum}}(t) \ge p_{\min}\delta$, therefore the condition $\delta \gtrsim \frac{\log n}{p_{\min}}$ implies $p_{\delta,\text{sum}}(t) \gtrsim \log n$, as well as $p_\delta(t) \gtrsim 1$ (due to Proposition 1), thus satisfying the requirements of Theorem 2. Additionally, $\delta$ is required to satisfy $\delta \le T$, and also condition (3.4), i.e.,

$$ \sqrt{\frac{\log n}{Ln\delta p_{\min}}} + \frac{M\delta n}{T} \lesssim \frac{1}{b^{\frac{7}{2}}(t)}. \qquad (4.10) $$

If $\delta$ satisfies the three aforementioned conditions, and if $n \gtrsim \log n$, we have with probability at least $1 - O(n^{-10})$ the $\ell_2$ bound

$$ \frac{\|\hat{\pi}(t) - \pi^*(t)\|_2}{\|\pi^*(t)\|_2} \lesssim \frac{M\delta nb^{\frac{7}{2}}(t)}{T} + \frac{b^{\frac{9}{2}}(t)}{\sqrt{Ln\delta p_{\min}}}. \qquad (4.11) $$

The optimal choice of $\delta \in (0, T]$ that minimizes the RHS of (4.11) is easily verified to be

$$ \delta^* = \min \left\{ \frac{(b(t)T)^{\frac{2}{3}}}{(2M)^{\frac{2}{3}}n(Lp_{\min})^{\frac{1}{3}}}, T \right\}. $$

Now it remains to ensure that $\delta^*$ satisfies the previously stated conditions on $\delta$. Clearly $\delta^* \le T$, and condition (4.10) are equivalent to the stated conditions on $T$ in the corollary. Hence provided $\delta \asymp \delta^*$, $T$ satisfies the stated conditions, and $n \gtrsim \log n$, we arrive at the stated $\ell_2$ bound in the corollary.
## 5 $\ell_\infty$-analysis of the spectral estimator

The main goal of this section is to present the steps of the proof of Theorem 3. We will follow the steps of the analysis carried out by Chen et al. [7], and adapt it to our setting. The proofs of all results from this section are provided in Appendix D.

As $\hat{\pi}(t)$ and $\pi^*(t)$ are stationary distributions of the transition matrices $\hat{P}(t)$ and $\bar{P}(t)$, then denoting $P_{.m}$ to be the $m^{th}$ column of a matrix $P$, it holds for all $m \in [n]$ that,

$$ 
\begin{aligned}
\hat{\pi}_m(t) - \pi_m^*(t) &= (\hat{\pi}(t)^\top \hat{P}(t))_m - (\pi^*(t)^\top \bar{P}(t))_m \\
&= \hat{\pi}(t)^\top \hat{P}_{.m}(t) - \pi^*(t)^\top \bar{P}_{.m}(t) \\
&= \pi^*(t)^\top (\hat{P}_{.m}(t) - \bar{P}_{.m}(t)) + (\hat{\pi}(t) - \pi^*(t))^\top \hat{P}_{.m}(t) \\
&= \underbrace{\pi^*(t)^\top (\hat{P}_{.m}(t) - \hat{P}_{.m}^*)(t)}_{:=I_0^m} + \underbrace{\pi^*(t)^\top (\hat{P}_{.m}^*(t) - \bar{P}_{.m}(t))}_{:=I_1^m} + \underbrace{(\hat{\pi}_m(t) - \pi_m^*(t)) \hat{P}_{mm}(t)}_{:=I_2^m} \\
&\quad + \sum_{j:j \neq m} (\hat{\pi}_j(t) - \pi_j^*(t)) \hat{P}_{jm}(t).
\end{aligned}
$$

Let us first focus on bounding $|I_0^m|$, $|I_1^m|$ and $|I_2^m|$. The last term will be treated carefully due to the occurrence of some statistical dependencies therein. Since $\hat{P}^*(t) = \mathbb{E}[\hat{P}(t)]$, hence $|I_0^m|$ can be bounded using Hoeffding's inequality.

**Lemma 6.** Suppose that $np_\delta(t) \ge c_0 \log n$ and $p_{\delta,\text{sum}}(t) \ge c_2 \log n$ for constants $c_0, c_2 \ge 1$ from Lemma 11. Then there exists a constant $C_1 \ge 1$ such that with probability at least $1 - O(n^{-9})$,

$$ |I_0^m| \le C_1 \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} \|\pi^*(t)\|_{\infty} \quad \forall m \in [n]. \quad (5.1) $$

The second term $I_1^m$ can be bounded easily using Assumption 1.

**Lemma 7.** With probability at least $1 - O(n^{-10})$, we have

$$ |I_1^m| \le \frac{2M\delta}{T} \|\pi^*(t)\|_{\infty} \quad \forall m \in [n]. \quad (5.2) $$

The next lemma shows that with high probability, $|I_2^m| \lesssim \|\hat{\pi}(t) - \pi^*(t)\|_{\infty}$ for all $m \in [n]$.

**Lemma 8.** Recall that $b_{\text{max},\delta}(t) = \max_{t \in \mathcal{N}_\delta(t)} b(t)$ and suppose that $np_\delta(t) \ge c_0 \log n$ and $p_{\delta,\text{sum}}(t) \ge c_2 \log n$ for constants $c_0, c_2 \ge 1$ from Lemma 11. Then there exists a constant $C_2 \ge 1$ such that it holds with probability at least $1 - O(n^{-9})$ that

$$ |I_2^m| \le \left( 1 - \frac{1}{12b_{\text{max},\delta}(t)} + C_2 \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} \right) \|\hat{\pi}(t) - \pi^*(t)\|_{\infty} \quad \forall m \in [n]. \quad (5.3) $$

The last term to bound is more difficult to handle due to the statistical dependency between $\hat{\pi}(t)$ and $\hat{P}(t)$. The idea is then to use the same “leave-one-out” trick as in [7] and introduce a new matrix $\hat{P}^{(m)}(t)$ with entries given by (for all $i \neq j$)

$$ \hat{P}_{ij}^{(m)}(t) = \begin{cases} \hat{P}_{ij}(t) & \text{if } i \neq m, j \neq m \\ \frac{p_{\delta}(t)}{d_{\delta}(t)} \frac{w_{t,j}^*}{w_{t,j}^* + w_{t,i}^*} & \text{if } i = m \text{ or } j = m. \end{cases} $$
To ensure that $\hat{P}^{(m)}(t)$ is a transition matrix, its diagonal entries are defined as

$$ \hat{P}_{ii}^{(m)}(t) = 1 - \sum_{k \neq i} \hat{P}_{ik}^{(m)}(t). $$

Here, the $m^{th}$ line and column of $\hat{P}(t)$ have been replaced by their expectation, unconditionally of the union graph $G_{\delta}(t)$. Let us denote $\hat{\pi}^{(m)}(t)$ to be the leading left eigenvector of $\hat{P}^{(m)}(t)$. This vector is now statistically independent of the connectivity and the pairwise comparison outputs involving the $m^{th}$ item, and one can reasonably expect that it is close to $\pi^{*}(t)$. Hence, we can decompose the last term as

$$ \sum_{j: j \neq m} (\hat{\pi}_j(t) - \pi_j^*(t)) \hat{P}_{jm}(t) = \underbrace{\sum_{j: j \neq m} (\hat{\pi}_j(t) - \hat{\pi}_j^{(m)}(t)) \hat{P}_{jm}(t)}_{:=I_3^m} + \underbrace{\sum_{j: j \neq m} (\hat{\pi}_j^{(m)}(t) - \pi_j^*(t)) \hat{P}_{jm}(t)}_{:=I_4^m} $$

To bound $|I_3^m|$, note that the Cauchy-Schwarz inequality implies

$$ |I_3^m| \leq \sqrt{\sum_{j: j \neq m} \hat{P}_{jm}(t)^2 \|\hat{\pi}(t) - \hat{\pi}^{(m)}(t)\|_2} $$

Since for all $j \neq m$, $\hat{P}_{jm}(t) \leq \frac{1}{d_{\delta}(t)}$, hence $|I_3^m| \leq \frac{\sqrt{d_{\delta, \max}(t)}}{d_{\delta}(t)} \|\hat{\pi}(t) - \hat{\pi}^{(m)}(t)\|_2$. The important step now is to bound $\|\hat{\pi}(t) - \hat{\pi}^{(m)}(t)\|_2$ since $d_{\delta, \max}(t) \leq d_{\delta}(t)$ w.h.p (when $G_{t'} \sim G(n, p(t'))$ for all $t' \in \mathcal{N}_{\delta}(t)$) due to Lemma 11. This is shown in the following lemma.

**Lemma 9.** Suppose that $np_{\delta}(t) \geq c_0 \log n$ and $p_{\delta, \text{sum}}(t) \geq c_2 \log n$ for constants $c_0, c_2 \geq 1$ from Lemma 11. Then there exist constants $C_3, C_6 \geq 1$ such that if

$$ 96b^{\frac{5}{2}}(t) \left( \frac{4M\delta}{T} + C_3 \sqrt{\frac{\log n}{np_{\delta}(t)}} \right) \leq \frac{1}{2}, \qquad (5.4) $$

then it holds with probability at least $1 - O(n^{-9})$ that for all $m \in [n]$,

$$ \|\hat{\pi}^{(m)}(t) - \hat{\pi}(t)\|_2 \leq 192b^{\frac{5}{2}}(t) \left( \frac{4M\delta}{T} + C_3 \sqrt{\frac{\log n}{Lnp_{\delta}(t)p_{\delta, \text{sum}}(t)}} \right) \|\pi^*(t)\|_{\infty} + \|\hat{\pi}(t) - \pi^*(t)\|_{\infty}. $$

Consequently, we have with probability at least $1 - O(n^{-9})$ that for all $m \in [n]$,

$$ |I_3^m| \leq 192 \frac{b^{\frac{5}{2}}(t)}{\sqrt{3np_{\delta}(t)}} \left( \frac{4M\delta}{T} + C_3 \sqrt{\frac{\log n}{Lnp_{\delta}(t)p_{\delta, \text{sum}}(t)}} \right) \|\pi^*(t)\|_{\infty} + \frac{1}{\sqrt{3np_{\delta}(t)}} \|\hat{\pi}(t) - \pi^*(t)\|_{\infty}. \qquad (5.5) $$

Finally, we can bound $|I_4^m|$ using the statistical independence between $\hat{\pi}^{(m)}(t)$ and $\hat{P}_{.m}(t), \hat{P}_{.m}(t)$.

**Lemma 10.** Suppose that $np_{\delta}(t) \geq c_0 \log n$ and $p_{\delta, \text{sum}}(t) \geq c_2 \log n$ for constants $c_0, c_2 \geq 1$ from Lemma 11; $n \geq c_1 \log n$ for the constant $c_1 > 0$ in Theorem 1, and that condition (3.4) of Theorem 2 holds. Then there exist constants $C_7, C_8, C_9 \geq 1$ such that with probability at least $1 - O(n^{-9})$, we have for all $m \in [n]$,

$$ |I_4^m| \leq \left( C_7 \frac{Mn\delta b^{\frac{7}{2}}(t)}{T} + C_8 \frac{b^{\frac{5}{2}}(t) \max\left\{b^2(t), \frac{\log n}{\sqrt{np_{\delta}(t)p_{\delta, \text{sum}}(t)}}\right\}}{\sqrt{Lnp_{\delta}(t)p_{\delta, \text{sum}}(t)}} \right) \|\pi^*(t)\|_{\infty} + C_9 \sqrt{\frac{\log n}{np_{\delta}(t)}} \|\hat{\pi}(t) - \pi^*(t)\|_{\infty}. $$
## 5.1 Proof of Theorem 3

As seen in Section 5, the bound in Theorem 3 follows from the combination of the bounds on $I_0^m, I_1^m, I_2^m, I_3^m$, and $I_4^m$. Note that we can identify two types of terms in these bounds – those depending on $\|\pi^*(t)\|_{\infty}$ and the ones which depend on $\|\hat{\pi}(t) - \pi^*(t)\|_{\infty}$. Then, our bound can be written as

$$ \|\hat{\pi}(t) - \pi^*(t)\|_{\infty} \le \alpha \|\pi^*(t)\|_{\infty} + \beta \|\hat{\pi}(t) - \pi^*(t)\|_{\infty} $$

where $\beta < 1$, which in turn implies $\frac{\|\hat{\pi}(t) - \pi^*(t)\|_{\infty}}{\|\pi^*(t)\|_{\infty}} \le \frac{\alpha}{1-\beta}$.

We can lower bound the term $1 - \beta$ as

$$ \begin{aligned} 1 - \beta &= 1 - \left( 1 - \frac{1}{12b_{\max,\delta}(t)} + C_2 \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} + \frac{1}{\sqrt{np_\delta(t)}} + C_9 \sqrt{\frac{\log n}{np_\delta(t)}} \right) \\ &= \frac{1}{12b_{\max,\delta}(t)} - C_2 \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} - \frac{1}{\sqrt{np_\delta(t)}} - C_9 \sqrt{\frac{\log n}{np_\delta(t)}} \\ &\ge \frac{1}{12b_{\max,\delta}(t)} - \tilde{C} \sqrt{\frac{\log n}{np_\delta(t)}} \quad (\text{since } p_{\delta,\text{sum}}(t) \ge \log n) \end{aligned} \quad (5.6) $$

for some constant $\tilde{C} \ge 1$.

Concerning the terms in $\alpha$, one can divide them in two groups depending on whether they depend on $T$ or not. The sum of the terms depending on $T$ is

$$ \frac{2M\delta}{T} + 768 \frac{M\delta b^{\frac{5}{2}}(t)}{T\sqrt{3np_\delta(t)}} + C_7 \frac{Mn\delta b^{\frac{7}{2}}(t)}{T} \le \tilde{C}_6 \frac{Mn\delta b^{\frac{7}{2}}(t)}{T} $$

for some constant $\tilde{C}_6 \ge 1$. Furthermore, the sum of the terms which are independent of $T$ is

$$ \begin{aligned} & C_1 \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} + 192C_3 b^{\frac{5}{2}}(t) \sqrt{\frac{\log n}{3Lp_{\delta,\text{sum}}(t)n^2 p_\delta^2(t)}} + C_8 \frac{b^{\frac{5}{2}}(t) \max\left\{b^2(t), \frac{\log n}{\sqrt{np_\delta(t)}}\right\}}{\sqrt{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} \\ & \le \left( C_1 + 192C_3 \frac{b^{\frac{5}{2}}(t)}{\sqrt{3np_\delta(t)}} + C_8 \frac{b^{\frac{5}{2}}(t)}{\sqrt{\log n}} \max\left\{b^2(t), \frac{\log n}{\sqrt{np_\delta(t)}}\right\} \right) \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} \\ & \le \tilde{C}_5 \underbrace{\left( 1 + \frac{b^{\frac{5}{2}}(t)}{\sqrt{\log n}} \max\left\{b^2(t), \frac{\log n}{\sqrt{np_\delta(t)}}\right\} \right)}_{=: \gamma_{\delta,n}(t)} \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} \quad (\text{since } np_\delta(t) \ge \log n) \end{aligned} $$

for some constant $\tilde{C}_5 \ge 1$. Thus we arrive at the bound

$$ \alpha \le \tilde{C}_5 \gamma_{n,\delta}(t) \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} + \tilde{C}_6 \frac{Mn\delta b^{\frac{7}{2}}(t)}{T}. \quad (5.7) $$

Combining (5.6) and (5.7), we readily arrive at the stated bound in Theorem 3.
## 5.2 Proof of Corollary 2

As in the proof of Corollary 1, since $p_{\delta,\text{sum}}(t) \ge p_{\min}\delta$, the condition $\delta \gtrsim \frac{\log n}{p_{\min}}$ implies $p_{\delta,\text{sum}}(t) \gtrsim \log n$ and $p_{\delta}(t) \gtrsim 1$. Moreover, $\delta$ has to satisfy $\delta \le T$ and condition (3.6), i.e.

$$ \frac{M\delta}{T} + \sqrt{\frac{\log n}{n}} \lesssim \frac{1}{b^{\frac{5}{2}}(t)}. $$

If $\delta$ satisfies these assumptions and if $n \gtrsim \log n$, we have with probability at least $1 - O(n^{-9})$,

$$ \frac{\|\hat{\pi}(t) - \pi^*(t)\|_{\infty}}{\|\pi^*(t)\|_{\infty}} \lesssim \frac{b_{\max}}{1 - b_{\max}\sqrt{\frac{\log n}{n}}} \left( \gamma_n(t) \sqrt{\frac{\log n}{Ln\delta p_{\min}}} + \frac{Mn\delta b^{\frac{7}{2}}(t)}{T} \right). \quad (5.8) $$

The optimal choice of $\delta \in (0, T]$ minimizing the RHS of (5.8) is given by

$$ \delta^* = \min \left\{ \frac{(\gamma_n(t)T)^{\frac{2}{3}}(\log n)^{\frac{1}{3}}}{(2M)^{\frac{2}{3}}nb^{\frac{7}{3}}(t)(Lp_{\min})^{\frac{1}{3}}}, T \right\}. $$

It now remains to check that $\delta^*$ satisfies the aforementioned conditions on $\delta$. Clearly $\delta^* \le T$ while $\delta^* \gtrsim \frac{\log n}{p_{\min}}$ and condition (3.6) are ensured for the stated conditions on $T$. Hence, for $\delta \asymp \delta^*$, if $T$ satisfies the stated conditions and $n \gtrsim \log n$, we arrive at the stated expression for the $\ell_{\infty}$ error bound in the corollary.

# 6 Experiments

We now empirically evaluate the performance of our method via numerical experiments<fnref n="1" /> on synthetic data, and on a real dataset. We will in particular compare our method with the MLE approach of Bong et al. [2], as it is for now the only other method we are aware of that theoretically analyzes a 'smoothly evolving' dynamic BTL model as us.

## 6.1 Synthetic data

The synthetic data is generated as follows.

1. For $i \in [n]$, we simulate the strength of item $i$ across the grid $w_i^* \in \mathbb{R}^{|T|}$ using the same gaussian process $GP(\mu_i, \Sigma_i)$ as in [2].

* $\mu_i = (\mu_i(0), \dots, \mu_i(T))$ with $\mu_i(t) \sim N(0, 0.1)$ for all $t \in T$, and

* $\Sigma_i$ is a Toeplitz symmetric matrix, with its coefficients in the first row defined as $\Sigma_{i,1t} = 1 - \frac{t}{T+1}$.

We then define $w_i^* = \exp(GP(\mu_i, \Sigma_i)) \in \mathbb{R}^{|T|}$ for each $i \in [n]$.

2. For all $t \in T$, we simulate an Erdös-Renyi comparison graph $G(n, p(t))$ with $p(t)$ chosen randomly from the interval $[\frac{1}{n}, \frac{\log n}{n}]$. We check that the union graph of all the data on the grid $T$ is connected. Indeed, it is a sufficient condition for the existence at all times $t \in [0, 1]$ of a $\delta$ such that the union graph $G_{\delta}(t)$ is connected (which is required for the ranking recovery).

<footnote n="1">¹Code available at: https://github.com/karle-eglantine/Dynamic_Rank_Centrality</footnote>
3. For all $t \in \mathcal{T}$, for all $1 \le i < j \le n$ and for all $l \in [L]$, we draw the outcomes of the comparisons as $y_{ij}^{(l)}(t) \sim \mathcal{B}\left(\frac{w_{t,i}^*}{w_{t,i}^*+w_{t,j}^*}\right)$ and define $y_{ji}^{(l)}(t) = 1 - y_{ij}^{(l)}(t)$.

Starting with an initial value of $\delta$ as in Corollary 1 ($\delta \approx T^{\frac{2}{3}}$), we increase its value till the union graph $G_{\delta}(t)$ is connected. We then recover the weight vectors $w_t^*$ for all $t \in \mathcal{T}$ using Algorithm 1. This process is repeated over 30 Monte Carlo runs.

Apart from the MLE approach [2], we also evaluate against an adaptation of the simple Borda Count method from the static setting to the dynamic setup.

**Borda Count.** This method, analysed by Ammar and Shah [1] in the static case, gives a score to each item based on its win rate. To estimate the scores at time $t$ in our dynamic setup, we compute the win rate of each item $i$ using the neighborhood $\mathcal{N}_{\delta}(t)$ as

$$s(t, i) = \frac{\sum_{t' \in \mathcal{N}_{\delta}(t)} \text{Number of times } i \text{ has won at time } t'}{\sum_{t' \in \mathcal{N}_{\delta}(t)} \text{Number of times } i \text{ has been compared at time } t'}$$

The scores $(s(t, i))_{i \in [n]}$ then yield a ranking of the items at time $t$ (the higher the winrate for an item, the better its rank).

**Ranking error.** In order to compare the rankings produced by these three methods, we compute for all them an estimation error with the error metric defined by Negahban [23]. Let $\pi^* = \frac{w^*}{\|w^*\|_1}$ denote the normalized true weight vector and $\sigma$ denote an estimated ranking, with $\sigma_i < \sigma_j \Leftrightarrow i$ is better than $j$. The error metric is then defined<fnref n="2" /> as follows.

$$D_{\pi^*}(\sigma) = \sqrt{\frac{1}{2n\|\pi^*\|_2^2} \sum_{i<j} (\pi_i^* - \pi_j^*)^2 \mathbf{1}_{(\pi_i^* - \pi_j^*)(\sigma_i - \sigma_j) > 0}.$$

It has been shown in [23, Lemma 1] that this error criterion is related to the $\ell_2$ error that we have bounded in Theorem 1. If the ranking $\sigma$ comes from a weight vector $\hat{\pi}$, then

$$D_{\pi^*}(\sigma) \le \frac{\|\pi^* - \hat{\pi}\|_2}{\|\pi^*\|_2}. \qquad (6.1)$$

Although (6.1) doesn't necessarily require $\hat{\pi}$ to satisfy $\|\hat{\pi}\|_1 = 1$, we will impose this in what follows since $\hat{\pi}$ will be the strength estimates returned by the methods being compared.

**Interpreting the results.** In Figure 1, we consider $T$ ranging from 10 to 150, and fix $n = 100$, $L = 5$. Figure 1a shows the evolution of the mean ranking error $D_{\pi^*}(\sigma)$ as a function of $T$, where for each $T$ the mean is taken across all time instants in $\mathcal{T}$, and all Monte Carlo runs. This is shown for our Algorithm 1 (dubbed DRC for Dynamic Rank Centrality), as well as the MLE and Borda Count methods, with $\sigma, \pi$ denoting the estimated ranks and weights by these algorithms. In Figure 1b, we plot the mean $\ell_2$ error $\frac{\|\pi^* - \hat{\pi}\|_2}{\|\pi^*\|_2}$ versus $T$ for only DRC and the MLE since Borda Count is not designed for recovering the latent weight vector $w^*$. One can observe that $\ell_2$ error decreases with $T$ for both DRC and MLE, which is consistent with the theoretical results for these methods developed in the present paper, and in [2]. Using (6.1), this implies that $D_{\pi^*}(\sigma)$ should decrease with $T$ which is what we observe in Figure 1a. One can note that both MLE and DRC have similar

<footnote n="2">²The quantities $\pi^*$, $\sigma$ and $\hat{\pi}$ are obviously defined at the time instant $t$ where the estimation is being carried out, however we suppress the dependence on $t$ for ease of exposition.</footnote>
performance for both error criteria. The Borda Count method performs well for rank recovery, as its error is only slightly worse than the other methods for $n = 100$. The error bars in Figure 1b illustrate that the variance of the errors typically decreases with $T$. In Figure 2, we show the evolution of the errors for $n = 400$ with the other parameters chosen as in Figure 1. We note that the observations made for the case $n = 100$ still hold and the error curves are lower than those for $n = 100$. Moreover, the variance of the errors in Figure 2b are slightly smaller than in Figure 1b.

<figure><img src="image_3.png" /><figcaption>**Figure 1:** Evolution of estimation errors with $T$ for Dynamic Rank Centrality, the MLE and Borda Count method for $n = 100$. The results are averaged over the grid $T$ as well as 30 Monte Carlo runs.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Figure 2:** Evolution of estimation errors with $T$ for Dynamic Rank Centrality, the MLE and Borda Count method for $n = 400$. The results are averaged over the grid $T$ as well as 30 Monte Carlo runs.</figcaption></figure>
## 6.2 Real dataset: NFL data

We now evaluate our method on a real dataset which consists of the results of National Football League (NFL) games for each season between 2009 and 2015, that are available in the `nflWAR` package [30]. The aim is to recover the ranking of the $n = 32$ teams at the end of each season, which contains $T = 16$ rounds. The dataset is hence composed of 16 comparison graphs with 32 nodes each, and comparison outcomes $(y_{ij}(t))_{i,j,t}$, where $y_{ij}(t) = 1$ if team $j$ beat team $i$ during in round $t$. We fit our model to this data by estimating the underlying strengths (and therefore the ranks) of the teams, at the end of a season. To do so, we tune the parameter $\delta$ using a Leave-One-Out Cross-Validation (LOOCV) procedure, described below.

1. Fix a list of potential values of $\delta$ and compute for each of them the associated estimates of the strength $\hat{\pi}$ at the end of the season.

2. For every possible values of $\delta$, repeat the following steps several times.

* Select randomly a game during the season, identified by a time $t$ and the pair of compared teams $\{i, j\}$.

* Consider the dataset obtained by removing the outcome $y_{ij}(t)$ of this game. Use this dataset to compute the estimated strength $\tilde{\pi}(t)$ at this time $t$.

* Compute the prediction error $\|y_{ij}(t) - \frac{\tilde{\pi}_j(t)}{\tilde{\pi}_j(t) + \tilde{\pi}_i(t)}\|_2$.

We then compute the mean of these prediction errors for each value of $\delta$.

3. Finally, we select $\delta^*$ which has the smallest mean prediction error, and take as our estimate of the strengths at the end of the season the associated vector $\hat{\pi}$.

An analogous LOOCV procedure is performed in the MLE approach [2]. We then compare our estimated ranking with the ones obtained by the MLE method, the Borda Count method and also with the ELO ratings. The latter are reputed to be relevant estimations of the teams qualities and are also openly available [22]. Figure 3 contains the estimation of the top 10 teams by all these methods, for the seasons 2011 to 2015. Considering the ELO ranks as the true ranks, one can observe that the DRC and MLE methods perform similarly, as a majority of the top 10 teams are recovered for each season. The Borda Count does not perform as well as it did on synthetic data; We can observe that it recovers the same ranks for several teams. This is explained by the fact that in this dataset, each team plays a small number of games, and thus the win rates take a finite (and small) number of values. However, it still recovers a large fraction of the top 10 teams.

# 7 Discussion

We now provide a detailed discussion with closely related work for dynamic ranking, and conclude with future directions for research.

## 7.1 Related work

As mentioned in Section 1, existing theoretical results for the dynamic ranking setup are limited, and the only works we are aware of are the recent results of Bong et al. [2] and of Li and Wakin [15]. We now discuss these two results in more detail.
<table><thead><tr><th rowspan="2">Ranks</th><th colspan="4">2011</th><th colspan="4">2012</th></tr><tr><th>ELO</th><th>DRC</th><th>MLE</th><th>Borda</th><th>ELO</th><th>DRC</th><th>MLE</th><th>Borda</th></tr></thead><tbody><tr><td>1</td><td>GB</td><td>GB</td><td>GB</td><td>GB</td><td>NE</td><td>NE</td><td>HOU</td><td>DEN</td></tr><tr><td>2</td><td>NE</td><td>BAL</td><td>DET</td><td>SF,NE,NO</td><td>DEN</td><td>SF</td><td>ATL</td><td>NE</td></tr><tr><td>3</td><td>NO</td><td>SF</td><td>NO</td><td></td><td>GB</td><td>DEN</td><td>SF</td><td>IND,ATL,GB</td></tr><tr><td>4</td><td>PIT</td><td>NE</td><td>NE</td><td></td><td>SF</td><td>SEA</td><td>CHI</td><td></td></tr><tr><td>5</td><td>BAL</td><td>NO</td><td>BUF</td><td>BAL,PIT</td><td>ATL</td><td>GB</td><td>GB</td><td></td></tr><tr><td>6</td><td>SF</td><td>PIT</td><td>SF</td><td></td><td>SEA</td><td>ATL</td><td>NE</td><td>SEA,HOU,WAS,SF</td></tr><tr><td>7</td><td>ATL</td><td>KC</td><td>CHI</td><td>DET,ATL,HOU</td><td>NYG</td><td>IND</td><td>DEN</td><td></td></tr><tr><td>8</td><td>PHI</td><td>NYG</td><td>TB</td><td></td><td>CIN</td><td>NYG</td><td>SEA</td><td></td></tr><tr><td>9</td><td>SD</td><td>HOU</td><td>OAK</td><td></td><td>BAL</td><td>STL</td><td>BAL</td><td></td></tr><tr><td>10</td><td>HOU</td><td>DET</td><td>DAL</td><td>TEN,NYG,CIN</td><td>HOU</td><td>WAS</td><td>IND</td><td>BAL,CHI,MIN,CIN</td></tr></tbody></table>

<table><caption>**Figure 3:** The top 10 teams for seasons 2011 to 2015, using ELO ranks, DRC, the MLE, and Borda Count. Teams highlighted in yellow for a particular recovery method are teams appearing in the top 10 list for ELO rankings. Teams highlighted in green are recovered at the same rank as in the ELO rankings.</caption><thead><tr><th rowspan="2">Ranks</th><th colspan="4">2013</th><th colspan="4">2014</th><th colspan="4">2015</th></tr><tr><th>ELO</th><th>DRC</th><th>MLE</th><th>Borda</th><th>ELO</th><th>DRC</th><th>MLE</th><th>Borda</th><th>ELO</th><th>DRC</th><th>MLE</th><th>Borda</th></tr></thead><tbody><tr><td>1</td><td>SEA</td><td>CAR</td><td>SEA</td><td>CAR</td><td>SEA</td><td>SEA</td><td>DEN</td><td>SEA,NE,DEN,GB,DAL</td><td>SEA</td><td>CAR</td><td>DAR</td><td>CAR</td></tr><tr><td>2</td><td>SF</td><td>SEA</td><td>DEN</td><td>SEA,SF</td><td>NE</td><td>DEN</td><td>ARI</td><td></td><td>CAR</td><td>ARI</td><td>DEN</td><td>ARI</td></tr><tr><td>3</td><td>NE</td><td>NO</td><td>NO</td><td></td><td>DEN</td><td>NE</td><td>NE</td><td></td><td>ARI</td><td>DEN</td><td>NE</td><td>NE,DEN,CIN</td></tr><tr><td>4</td><td>DEN</td><td>SF</td><td>KC</td><td>ARI,PIT,CIN,PHI,DEN,NE,IND</td><td>GB</td><td>ARI</td><td>SEA</td><td></td><td>KC</td><td>CIN</td><td>CIN</td><td></td></tr><tr><td>5</td><td>CAR</td><td>ARI</td><td>SF</td><td></td><td>DAL</td><td>DAL</td><td>DAL</td><td></td><td>DEN</td><td>PIT</td><td>ARI</td><td></td></tr><tr><td>6</td><td>CIN</td><td>NYJ</td><td>NE</td><td></td><td>PIT</td><td>GB</td><td>GB</td><td>ARI,PIT,DET,IND</td><td>NE</td><td>NE</td><td>GB</td><td>MIN,KC</td></tr><tr><td>7</td><td>NO</td><td>PHI</td><td>IND</td><td></td><td>BAL</td><td>KC</td><td>PHI</td><td></td><td>PIT</td><td>SEA</td><td>MIN</td><td></td></tr><tr><td>8</td><td>ARI</td><td>CIN</td><td>CAR</td><td></td><td>IND</td><td>PIT</td><td>SD</td><td></td><td>CIN</td><td>MIN</td><td>KC</td><td>GB,NYJ,SEA,PIT</td></tr><tr><td>9</td><td>IND</td><td>NE</td><td>ARI</td><td></td><td>ARI</td><td>CIN</td><td>DET</td><td></td><td>GB</td><td>KC</td><td>PIT</td><td></td></tr><tr><td>10</td><td>SD</td><td>STL</td><td>CIN</td><td></td><td>CIN</td><td>SD</td><td>KC</td><td>BAL,PHI,CIN</td><td>MIN</td><td>ATL</td><td>SEA</td><td></td></tr></tbody></table>

**Comparison to Bong et al.[2].** The dynamic BTL model proposed by Bong et al. [2] is closely related to the one we presented in Section 2. They consider the logit version of the BTL model

$$ \text{logit}(\mathbb{P}(i \text{ beats } j \text{ at time } t)) = \beta_i^*(t) - \beta_j^*(t) $$

where $\beta^*(t)$ represents the vector of scores at time $t$ with $w^*(t) = \exp(\beta^*(t))$. The main differences with our model are that (a) the grid $\mathcal{T}$ can be non-uniform, and (b) the number of comparisons $L_{ij}(t)$ made for each pair $\{i, j\}$ at each time $t \in [0, 1]$ can vary. We assumed the grid $\mathcal{T}$ to be uniform only for simplicity, our analysis can be easily extended to handle the non-uniform setting as long as $\mathcal{T}$ is “sufficiently regular”. Similarly, one can easily extend our analysis to handle varying number of comparisons $L_{ij}(t)$; the assumption $L_{ij}(t) = L$ in our setup was only made for simplicity in the proofs as is typically done for the static setting (see for eg., [6, Remark 3]).

The pairwise comparison data at each time $t' \in \mathcal{T}$ are gathered into a matrix $X(t')$ where $X_{ij}(t')$ is the number of times $i$ beat $j$ at time $t'$. In order to use the temporal aspect of the data, they smooth the data using a kernel function; this is analogous to (weighted) averaging the data in a suitable neighborhood. More precisely, to estimate the ranks at time $t$, they first compute the smoothed data as

$$ \tilde{X}(t) = \sum_{t'=0}^{T} W_h(t', t)X(t') $$

where $W_h$ is a kernel function with bandwidth $h$. Then, $\beta^*$ is estimated by minimizing the negative
log-likelihood, i.e.,

$$ \underset{\beta: \sum_{i=1}^{n} \beta_i = 0}{\operatorname{argmin}} \quad \underset{\beta: \sum_{i=1}^{n} \beta_i = 0}{\operatorname{argmin}} \quad \left( \sum_{i,j: i \neq j} \frac{\tilde{X}_{ij}(t) \log(1 + \exp(\beta_j(t) - \beta_i(t)))}{\sum_{i,j: i \neq j} \tilde{X}_{ij}(t)} \right) \quad (7.1) $$

using a proximal gradient descent algorithm. If the matrix $|\tilde{X}(t)|$ is considered to be the adjacency matrix of a weighted directed graph, then the strong connectivity of this graph is sufficient for the unique existence of the solution of (7.1). The main assumptions needed for their theoretical analysis are the following.

* The probabilities $P(i \text{ beats } j \text{ at time } t)$ are Lipschitz functions of time $t \in [0, 1]$ for all $i \neq j \in [n]$ (same as Assumption 1).

* Each pair of teams $\{i, j\}$ has been compared at least at one time point $t' \in \mathcal{T}$. Translated to our notation, this means that the union graph $\cup_{t' \in \mathcal{T}} G_{t'}$ is complete.

We remark that this last assumption is a stronger assumption than the connectivity assumption on the union graph $G_{\delta}(t)$ that we required in our analysis.

Bong et al. derive bounds in the $\ell_{\infty}$ norm for estimating $\beta^*$. For a bandwidth $h \gtrsim (\log n/T)^{1/3}$, denoting $\hat{\beta}(t)$ to be the MLE estimator (i.e., the solution of (7.1)), it is shown [2, Theorem 5.2] that with high probability,

$$ \|\hat{\beta}(t) - \beta^*(t)\|_{\infty} \lesssim \delta_h(t) + \left(\frac{\log n}{nT}\right)^{\frac{1}{3}} $$

where $\delta_h(t)$ is denoted to be a discrepancy parameter in [2], and is small when all the teams play the same number of games at all times. However, the dependence of $\delta_h$ on $h$ is not made explicit. Note that they also recover the $T^{-1/3}$ rate for pointwise estimation of Lipschitz functions. Additionally, they also derive a rate for the uniform error. For $h \gtrsim (\log(nT^3)/T)^{1/3}$, it is shown [2, Theorem 5.3] that with high probability,

$$ \sup_{t \in [0,1]} \|\hat{\beta}(t) - \beta(t)\|_{\infty} \lesssim \sup_{t \in [0,1]} \delta_h(t) + \left(\frac{\log(nT)}{T}\right)^{\frac{1}{3}}. $$

While we provide pointwise estimation error bounds for any given $t \in [0, 1]$, it is also possible to extend our results to obtain error bounds holding uniformly over all $t \in [0, 1]$. The first main idea here would be to observe that there are $O(T)$ different number of neighborhoods in the construction of $\hat{P}(t)$ in (2.7), which implies that there are $O(T)$ different possible values of $\hat{\pi}(t)$ over all $t \in [0, 1]$. Secondly, one can verify that $\pi^*(t)$ is a Lipschitz function of $t$. These two facts together with a union bound argument can be used to establish $\ell_2$ and $\ell_{\infty}$ bounds holding uniformly over $t \in [0, 1]$, with probability at least $1 - O(Tn^{-c})$ where $c$ is a suitably large constant. In our analysis, we had taken $c$ to be 10 (resp. 9) for the $\ell_2$ (resp. $\ell_{\infty}$) analysis, but it can be any suitably large value, at the expense of worsening the other constants in the accompanying theorems in Section 3.

**Comparison to Li and Wakin [15].** The model introduced by Li and Wakin [15] aims at recovering a pairwise comparison matrix $X(T)$ at a time $T$ from noisy linear measurements of the matrices $X(t)$ for $t \in [T]$. For each pair of items $\{i, j\}$, $X_{ij}(t)$ quantifies the preference of item $j$
over item i at time t with $X_{ij}(t) > 0$ if item j is preferred to item i. One can for example consider $X_{ij}(t)$ to be the goal difference during a football game between the teams j and i at time t. Note that these data contain more information than the simple binary outcomes we use in the dynamic BTL model. A specificity of this model is that the preferences over the collection of items are supposed to depend on $2r$ latent factors, captured in the matrices $S_t \in \mathbb{R}^{n \times r}$ and $Q \in \mathbb{R}^{n \times r}$. The model considered in [15] assumes that

$$X(t) = S_t Q^\top - Q S_t^\top \quad \forall t \in [T].$$

The matrix $Q$ contains information on factors that do not depend on time. The time-dependent data are contained in the matrix $S_t$ and evolve under the following generative model.

$$S_t = S_{t-1} + E_t \quad \forall t \in [T]$$

where $E_t$ is an innovation matrix with i.i.d. centered Gaussian entries. The intuition behind this model comes from the case of a single factor ($r = 1$) with $Q \in \mathbb{R}^n$ being the all ones vector. In this case, for each pair of items $\{i, j\}$,

$$X_{ij}(t) = s_{t,i} - s_{t,j} \quad \forall t \in [T].$$

The outcomes of the comparisons then, as in the BTL, only depend on the strength of each item at time $t$. Then by recovering the matrix $X(T)$, one can subsequently also derive a ranking of the items (see for eg., [8]).

Denoting $M$ to be the number of measurements available at each time $t$, they show a bound on the estimation error $\|X(T) - \hat{X}(T)\|_F$ that holds with high probability and decreases with $M$, where $\hat{X}(T)$ is computed as the solution of an optimization problem. It is however unclear how one can then derive a ranking of the items from this estimated comparison matrix $\hat{X}(T)$ for $r > 1$.

## 7.2 Future work

For future work it would be interesting to extend our study to other models. Indeed, the BTL model is classic for the ranking problem but suffers from some limitations.

1. *Independence.* One of the main assumptions of the dynamic BTL model we consider is the independence of the outcomes for all comparisons at a given time $t$, and the independence across different time points. The latter assumption can in particular be questioned, as the choices of one user across time are typically not going to be independent. It would be interesting to model these dependencies across time, for instance, by modelling them as a Markov process.

2. *Home effect.* Another model, used in the study of sports tournaments by Cattelan et al. [5] adapts the BTL to include a home effect, which is beneficial to the hosting team. Indeed, because of the familiar environment and the supporting public, a team is more likely to win if they play in their stadium. It would be interesting to theoretically analyze such a model.

3. *Non-parametric dynamic models.* Considering a parametric model as the BTL model has limitations as it assumes that the preference outcomes depend only on one parameter $w^*$, which can be seen as a strong transitivity constraint [27]. That is why one can instead use non-parametric models as described in [27, 28, 25], where only a certain transitivity assumption is made on the comparison matrix. In particular, the BTL model belongs to this class of models. It has been shown by Shah et al. [27] that in the static case, the matrix of
probabilities can be estimated at the same rate as in a parametric model. As such extensions
can be considered in the static setting, it would be interesting to adapt these models to the
dynamic case and study their performance.

References

[1] Ammar Ammar and Devavrat Shah. Ranking: Compare, don't score. In *2011 49th Annual Allerton Conference on Communication, Control, and Computing (Allerton)*, pages 776-783. IEEE, 2011.

[2] Heejong Bong, Wanshan Li, Shamindra Shrotriya, and Alessandro Rinaldo. Nonparametric estimation in the dynamic bradley-terry model. In *International Conference on Artificial Intelligence and Statistics*, pages 3317-3326. PMLR, 2020.

[3] Ralph Allan Bradley and Milton E Terry. Rank analysis of incomplete block designs: I. the method of paired comparisons. *Biometrika*, 39(3/4):324-345, 1952.

[4] Anne E Carpenter, Thouis R Jones, Michael R Lamprecht, Colin Clarke, In Han Kang, Ola Friman, David A Guertin, Joo Han Chang, Robert A Lindquist, Jason Moffat, et al. Cellprofiler: image analysis software for identifying and quantifying cell phenotypes. *Genome biology*, 7(10):1-11, 2006.

[5] Manuela Cattelan, Cristiano Varin, and David Firth. Dynamic bradley-terry modelling of sports tournaments. *Journal of the Royal Statistical Society: Series C (Applied Statistics)*, 62(1):135-150, 2013.

[6] Pinhan Chen, Chao Gao, and Anderson Y Zhang. Partial recovery for top-k ranking: Optimality of mle and sub-optimality of spectral method. *arXiv preprint arXiv:2006.16485*, 2020.

[7] Yuxin Chen, Jianqing Fan, Cong Ma, and Kaizheng Wang. Spectral method and regularized mle are both optimal for top-k ranking. *The Annals of Statistics*, 47(4):2204-2235, Aug 2019.

[8] Alexandre d'Aspremont, Mihai Cucuringu, and Hemant Tyagi. Ranking and synchronization from pairwise measurements via SVD. *J. Mach. Learn. Res.*, 22:19:1-19:63, 2021.

[9] Ludwig Fahrmeir and Gerhard Tutz. Dynamic stochastic models for time-dependent ordered paired comparison systems. *Journal of the American Statistical Association*, 89(428):1438-1449, 1994.

[10] Mark E Glickman and Hal S Stern. A state-space model for national football league scores. *Journal of the American Statistical Association*, 93(441):25-35, 1998.

[11] Pierre-Emmanuel Jabin and Stéphane Junca. A continuous model for ratings. *SIAM Journal on Applied Mathematics*, 75(2):420-442, 2015.

[12] Dietmar Jannach, Paul Resnick, Alexander Tuzhilin, and Markus Zanker. Recommender systems—beyond matrix completion. *Communications of the ACM*, 59(11):94-102, 2016.

[13] David A Levin and Yuval Peres. *Markov chains and mixing times*, volume 107. American Mathematical Soc., 2017.
[14] Lel Li and Karl Kim. Estimating driver crash risks based on the extended bradley-terry model: an induced exposure method. *Journal of the Royal Statistical Society: Series A (Statistics in Society)*, 163(2):227-240, 2000.

[15] Shuang Li and Michael B Wakin. Recovery guarantees for time-varying pairwise comparison matrices with non-transitivity. *arXiv preprint arXiv:2106.09151*, 2021.

[16] Michael J Lopez, Gregory J Matthews, and Benjamin S Baumer. How often does the best team win? a unified approach to understanding randomness in north american sport. *The Annals of Applied Statistics*, 12(4):2483-2516, 2018.

[17] R Duncan Luce. *Individual choice behavior: A theoretical analysis*. Courier Corporation, 2012.

[18] JNS Matthews and KP Morris. An application of bradley-terry-type models to the measurement of pain. *Journal of the Royal Statistical Society: Series C (Applied Statistics)*, 44(2):243-255, 1995.

[19] Lucas Maystre, Victor Kristof, and Matthias Grossglauser. Pairwise comparisons with flexible time-dynamics. In *Proceedings of the 25th ACM SIGKDD International Conference on Knowledge Discovery & Data Mining*, pages 1236-1246, 2019.

[20] Michael Mitzenmacher and Eli Upfal. *Probability and computing: Randomization and probabilistic techniques in algorithms and data analysis*. Cambridge university press, 2017.

[21] Shun Motegi and Naoki Masuda. A network-based dynamical ranking system for competitive sports. *Scientific reports*, 2(1):1-7, 2012.

[22] Paine N. Nfl elo ratings are back! https://fivethirtyeight.com/features/nfl-elo-ratings-are-back/, 2015.

[23] Sahand Negahban, Sewoong Oh, and Devavrat Shah. Rank centrality: Ranking from pairwise comparisons. *Operations Research*, 65(1):266-287, 2017.

[24] A. Nemirovski. Topics in non-parametric statistics. *Ecole d'Eté de Probabilités de Saint-Flour*, 28:85, 2000.

[25] Ashwin Pananjady, Cheng Mao, Vidya Muthukumar, Martin J. Wainwright, and Thomas A. Courtade. Worst-case versus average-case design for estimation from partial pairwise comparisons. *The Annals of Statistics*, 48(2):1072-1097, 2020.

[26] Mark Rudelson, Roman Vershynin, et al. Hanson-wright inequality and sub-gaussian concentration. *Electronic Communications in Probability*, 18, 2013.

[27] Nihar Shah, Sivaraman Balakrishnan, Aditya Guntuboyina, and Martin Wainwright. Stochastically transitive models for pairwise comparisons: Statistical and computational issues. In *International Conference on Machine Learning*, pages 11-20. PMLR, 2016.

[28] Nihar B Shah and Martin J Wainwright. Simple, robust and optimal ranking from pairwise comparisons. *The Journal of Machine Learning Research*, 18(1):7246-7283, 2017.

[29] Joel A. Tropp. An introduction to matrix concentration inequalities. *Foundations and Trends in Machine Learning*, 8(1-2):1-230, 2015.
[30] Ronald Yurko, Samuel Ventura, and Maksim Horowitz. nflwar: a reproducible method for offensive player evaluation in football. Journal of Quantitative Analysis in Sports, 15(3):163–183, 2019.

## A Technical tools

We collect here some technical results that are used for proving our main results. We begin by recalling the following result from [7].

**Theorem 4** ([7, Theorem 8]). Suppose that $P, \hat{P}, P^*$ are probability transition matrices with stationary distributions $\pi, \hat{\pi}, \pi^*$ respectively. We assume that $P^*$ represents a reversible Markov chain. When $\|P - \hat{P}\|_{\pi^*} < 1 - \max\{\lambda_2(P^*), -\lambda_n(P^*)\}$ holds, then we have

$$ \|\pi - \hat{\pi}\|_{\pi^*} \le \frac{\|\pi^\top(P - \hat{P})\|_{\pi^*}}{1 - \max\{\lambda_2(P^*), -\lambda_n(P^*)\} - \|P^* - \hat{P}\|_{\pi^*}}. \quad (A.1) $$

Next, we recall some useful concentration results from probability starting with the classical Chernoff bound.

**Theorem 5 (Chernoff bound, [20, Corollary 4.6]).** Let $X_1, \dots, X_n$ be independent Bernoulli variables with $\mathbb{P}(X_i = 1) = p_i$. Let $X = \sum_{i=1}^n X_i$ and $\mu = \mathbb{E}[X]$. For $\varepsilon \in (0, 1)$,

$$ \mathbb{P}(|X - \mu| \ge \varepsilon\mu) \le 2 \exp \left( -\frac{\mu\varepsilon^2}{3} \right). $$

**Theorem 6 (Hoeffding's inequality).** Let $X_1, \dots, X_n$ be a sequence of independent random variables where $X_i \in [a_i, b_i]$ for each $i \in [n]$, and $S_n = \sum_{i=1}^n X_i$. Then

$$ \mathbb{P}(|S_n - \mathbb{E}[S_n]| \ge t) \le 2 \exp \left( -\frac{2t^2}{\sum_{i=1}^n (b_i - a_i)^2} \right). $$

**Theorem 7 (Matrix Bernstein, [29]).** Let $Z_1, \dots, Z_n \in \mathbb{R}^{d_1 \times d_2}$ be independent, zero-mean random matrices, each satisfying (almost surely) $\|Z_i\|_2 \le B$. Then for any $s \ge 0$,

$$ \mathbb{P} \left( \left\| \sum_{i=1}^{n} Z_i \right\|_2 \ge s \right) \le (d_1 + d_2) \exp \left( -\frac{3s^2}{6\nu + 2Bs} \right) $$

where $\nu = \max \{ \|\mathbb{E}[\sum_{i=1}^n Z_i^\top Z_i]\|_2, \|\mathbb{E}[\sum_{i=1}^n Z_i Z_i^\top]\|_2 \}$.

## B Properties of Erdös-Renyi graphs

Recall that in our setup, we have $T + 1$ Erdös-Renyi graphs $G_{t'} \sim \mathcal{G}(n, p(t'))$ for $t' \in \mathcal{T}$. Moreover, for any given $t \in [0, 1]$, we have the union graph $G_{\delta}(t) \sim \mathcal{G}(n, p_{\delta}(t))$ with $p_{\delta}(t)$ as in (3.3). Also recall that $p_{\delta, \text{sum}}(t) := \sum_{t' \in \mathcal{N}_{\delta}(t)} p(t')$.

**Lemma 11.** Consider the events

$$ 1. \mathcal{A}_1 = \left\{ \frac{np_{\delta}(t)}{2} \le d_{\delta, \min}(t) \le d_{\delta, \max}(t) \le \frac{3np_{\delta}(t)}{2} \right\}, $$

$$ 2. \mathcal{A}_2 = \{\xi_{\delta}(t) > \frac{1}{2}\}, $$
$$3. \mathcal{A}_3 = \{|\mathcal{E}_\delta(t)| \le 2n^2p_\delta(t)\},$$

$$4. \mathcal{A}_4 = \left\{ \frac{p_{\delta, \text{sum}}(t)}{2} \le N_{\text{min}, \delta}(t) \le N_{\text{max}, \delta}(t) \le 2p_{\delta, \text{sum}}(t) \right\}.$$

Then the following is true.

1. There exists a constant $c_0 \ge 1$ such that if $p_\delta(t) \ge c_0 \frac{\log n}{n}$, then $\mathbb{P}(\mathcal{A}_i) \ge 1 - O(n^{-10})$ for $i = 1, 2, 3$.

2. There exists a constant $c_2 \ge 1$ such that if $p_{\delta, \text{sum}}(t) \ge c_2 \log n$, then $\mathbb{P}(\mathcal{A}_4) \ge 1 - O(n^{-10})$.

*Proof.* 1. This follows for $\mathcal{A}_1$ from the Chernoff bound (see Appendix A) with $\varepsilon = 1/2$, and the union bound. For $\mathcal{A}_2$, the statement follows directly by applying [23, Lemma 7] to $G_\delta(t)$. For $\mathcal{A}_3$, this is again a standard use of the Chernoff bound with $\mu = \binom{n}{2}p_\delta(t)$ and $\varepsilon = 1/2$.

2. Note that for any $i \ne j \in [n]$, $|\mathcal{N}_{ij,\delta}(t)|$ is a sum of independent Bernoulli random variables and $\mathbb{E}[|\mathcal{N}_{ij,\delta}(t)|] = \sum_{t' \in \mathcal{N}_{\delta}(t)} p(t') = p_{\delta, \text{sum}}(t)$. Then if $p_{\delta, \text{sum}}(t) \ge c_2 \log n$ for a suitably large constant $c_2 > 0$, it holds with probability at least $1 - O(n^{-12})$ that

$$\frac{p_{\delta, \text{sum}}(t)}{2} \le |\mathcal{N}_{ij, \delta}(t)| \le 2p_{\delta, \text{sum}}(t).$$

Now the union bound implies that

$$\mathbb{P}\left(\forall i \ne j : \frac{p_{\delta, \text{sum}}(t)}{2} \le |\mathcal{N}_{ij, \delta}(t)| \le 2p_{\delta, \text{sum}}(t)\right) \ge 1 - O(n^{-10}).$$

Since $N_{\text{min},\delta}(t) \ge \min_{i \ne j} |\mathcal{N}_{ij,\delta}(t)|$ and $N_{\text{max},\delta}(t) \le \max_{i \ne j} |\mathcal{N}_{ij,\delta}(t)|$ is always true, hence the statement follows. $\square$

**Proposition 1.** With $p_{\min}$ as in (3.5), we have for all $t \in [0, 1]$ that $p_\delta(t) \ge \min\{1, p_{\min}\delta\}$. Therefore $p_\delta(t) \ge \frac{\log n}{n}$ if $p_{\min}\delta \ge \frac{\log n}{n}$.

*Proof.* Starting with the definition of $p_\delta(t)$ in (3.3), we have

$$p_{\delta}(t) = 1 - \prod_{t' \in \mathcal{N}_{\delta}(t)} (1 - p(t')) \ge 1 - (1 - p_{\min})^{|\mathcal{N}_{\delta}(t)|} \ge 1 - e^{-p_{\min}|\mathcal{N}_{\delta}(t)|} \ge 1 - e^{-p_{\min}\delta} \quad (B.1)$$

where the last inequality follows uses the fact $|\mathcal{N}_\delta(t)| \ge \delta$ for all $t$. Now there are two cases to distinguish.

* If $p_{\min}\delta \le 1$, then $1 - e^{-p_{\min}\delta} \ge \frac{p_{\min}\delta}{2}$ using the standard inequality $1 - e^{-x} \ge x/2$ for $x \in [0, 1]$.

* If $p_{\min}\delta \ge 1$, then $1 - e^{-p_{\min}\delta} \ge 1 - e^{-1}$.

This concludes the proof. $\square$
# C Proofs of results in Section 4

## C.1 Proof of Lemma 1

Recall that the entries of $\Delta_1(t)$ are given by

$$ \Delta_{1,ij}(t) = \begin{cases} \frac{1}{d_{\delta}(t)|\mathcal{N}_{ij,\delta}(t)|} \sum_{t' \in \mathcal{N}_{ij,\delta}(t)} (y_{ij}^*(t') - y_{ij}^*(t)) & \text{if } \{i,j\} \in \mathcal{E}_{\delta}(t), \\ -\sum_{k \neq i} \Delta_{1,ik}(t) & \text{if } i = j, \\ 0 & \text{otherwise.} \end{cases} $$

Let us denote $D_1(t)$ to be the diagonal matrix containing the elements $\Delta_{1,ii}(t)$ and $\overline{\Delta}_1(t) = \Delta_1(t) - D_1(t)$. As $D_1(t)$ is diagonal, we then have

$$ \|\Delta_1(t)\|_2 \le \|D_1(t)\|_2 + \|\overline{\Delta}_1(t)\|_2 \le \max_i |\Delta_{1,ii}(t)| + \|\overline{\Delta}_1(t)\|_F. $$

Let us first bound $\|\overline{\Delta}_1(t)\|_F$. Assumption 1 implies that for any $\{i,j\} \in \mathcal{E}_{\delta}(t)$,

$$ |\overline{\Delta}_{1,ij}(t)| \le \frac{1}{d_{\delta}(t)|\mathcal{N}_{ij,\delta}(t)|} \sum_{t' \in \mathcal{N}_{ij,\delta}(t)} |y_{ij}^*(t') - y_{ij}^*(t)| \le \frac{1}{d_{\delta}(t)|\mathcal{N}_{ij,\delta}(t)|} \sum_{t' \in \mathcal{N}_{ij,\delta}(t)} M |t' - t| \le \frac{M\delta}{Td_{\delta}(t)}, $$

which in turn implies $\|\overline{\Delta}_1(t)\|_F \le 2 \frac{M|\mathcal{E}_{\delta}(t)|\delta}{Td_{\delta}(t)}$. In order to bound $\|D_1(t)\|_2$, we simply note that

$$ |D_{1,ii}(t)| = \left| -\sum_{j \neq i} \Delta_{1,ij}(t) \right| \le d_{\delta,\max}(t) \max_{j \neq i} |\Delta_{1,ij}(t)| \le d_{\delta,\max}(t) \frac{M\delta}{Td_{\delta}(t)} \le \frac{M\delta}{T} $$

since $d_{\delta}(t) \ge d_{\delta,\max}(t)$. Hence $\|D_1(t)\|_2 \le \frac{M\delta}{T}$, and so

$$ \|\Delta_1(t)\|_2 \le \frac{M\delta}{T} \left( 1 + 2 \frac{|\mathcal{E}_{\delta}(t)|}{d_{\delta}(t)} \right) \le \frac{M\delta}{T} \left( 1 + 2 \frac{|\mathcal{E}_{\delta}(t)|}{d_{\delta,\max}(t)} \right) \le 4 \frac{M\delta |\mathcal{E}_{\delta}(t)|}{Td_{\delta,\max}(t)} \quad (\text{since } 1 \le 2 \frac{|\mathcal{E}_{\delta}(t)|}{d_{\delta,\max}(t)}). $$

□

## C.2 Proof of Lemma 2

Our goal is to bound $\|\Delta(t)\|_2$ for the random matrix $\Delta(t) = \hat{P}(t) - \bar{P}(t)$, whose entries are

$$ \Delta_{ij}(t) = \begin{cases} \frac{1}{d_{\delta}(t)|\mathcal{N}_{ij,\delta}(t)|} \sum_{t' \in \mathcal{N}_{ij,\delta}(t)} (y_{ij}^*(t') - y_{ij}^*(t)) & \text{if } \{i,j\} \in \mathcal{E}_{\delta}(t), \\ -\sum_{k \neq i} \Delta_{ik}(t) & \text{if } i = j, \\ 0 & \text{otherwise.} \end{cases} $$

Once again $\Delta(t)$ can be decomposed as $\Delta(t) = D(t) + \overline{\Delta}(t)$ where $D(t)$ is a diagonal matrix and $\overline{\Delta}(t)$ contains its off-diagonal coefficients. Since

$$ \|\Delta(t)\|_2 \le \|D(t)\|_2 + \|\overline{\Delta}(t)\|_2, $$

we will bound $\|D(t)\|_2$ and $\|\overline{\Delta}(t)\|_2$ using ideas from the proof of [23, Lemma 3].
### C.2.1 Bound on $||D(t)||_2$

The matrix $D(t)$ is diagonal, hence $||D(t)||_2 = \max_{i \in [n]} |D_{ii}(t)|$ where the diagonal entries are given by

$$D_{ii}(t) = -\sum_{k \neq i} \sum_{t' \in \mathcal{N}_{ik,\delta}(t)} \sum_{l=1}^{L} \frac{1}{L |\mathcal{N}_{ik,\delta}(t)| d_{\delta}(t)} (y_{ik}^{(l)}(t') - y_{ik}^{*}(t')).$$

This is a sum of at most $Ld_{\delta,\max}(t)N_{\max,\delta}(t)$ independent and centered variables, which almost surely lie in the interval $[-\frac{1}{Ld_{\delta}(t)N_{\min,\delta}(t)}, \frac{1}{Ld_{\delta}(t)N_{\min,\delta}(t)}]$. Thus, by applying Hoeffding's inequality (see Theorem 6 in Appendix A) we obtain for any $s > 0$

$$\mathbb{P}(|D_{ii}(t)| > s) \le 2 \exp \left( - \frac{2s^2 L^2 d_{\delta}^2(t) N_{\min,\delta}^2(t)}{4Ld_{\delta,\max}(t) N_{\max,\delta}(t)} \right) = 2 \exp \left( - \frac{Ld_{\delta}^2(t) N_{\min,\delta}^2(t) s^2}{2N_{\max,\delta}(t) d_{\delta,\max}(t)} \right)$$

which implies via a union bound that

$$\mathbb{P}(\|D(t)\|_2 > s) \le 2n \exp \left( - \frac{Ld_{\delta}^2(t) N_{\min,\delta}^2(t) s^2}{2N_{\max,\delta}(t) d_{\delta,\max}(t)} \right).$$

Choosing $s = c_1 \sqrt{\frac{N_{\max,\delta}(t) d_{\delta,\max}(t) \log n}{Ld_{\delta}^2(t) N_{\min,\delta}^2(t)}}$ where $c_1 \ge 2$ is a constant, we have that

$$\|D(t)\|_2 \le c_1 \sqrt{\frac{N_{\max,\delta}(t) d_{\delta,\max}(t) \log n}{Ld_{\delta}^2(t) N_{\min,\delta}^2(t)}}$$

with probability at least $1 - 2n^{1-\frac{c_1^2}{2}}$.

### C.2.2 Bound on $||\bar{\Delta}(t)||_2$ when $d_{\delta,\max}(t)N_{\max,\delta}(t) \le \log n$

In order to bound $||\bar{\Delta}(t)||_2$, we first recall the standard inequality $||M||_2 \le \sqrt{||M||_1 ||M||_\infty}$, and also $||M||_1 = ||M^\top||_\infty$. Since $\bar{\Delta}(t)$ is skew-symmetric, it follows that

$$\|\bar{\Delta}(t)\|_2 \le \sqrt{\|\bar{\Delta}(t)\|_1 \|\bar{\Delta}(t)^\top\|_1} \le \|\bar{\Delta}(t)\|_1.$$

Our aim now is to bound $||\bar{\Delta}(t)||_1$. Recall that for every $i \neq j$,

$$\bar{\Delta}_{ij}(t) = \begin{cases} \frac{1}{d_{\delta}(t)|\mathcal{N}_{ij,\delta}(t)|} \sum_{t' \in \mathcal{N}_{ij,\delta}(t)} (y_{ij}(t') - y_{ij}^*(t')) & \text{if } \{i,j\} \in \mathcal{E}_{\delta}(t) \\ 0 & \text{otherwise.} \end{cases}$$

Now denoting for all $i$,

$$R_i(t) = \sum_{j \neq i} \sum_{t' \in \mathcal{N}_{ij,\delta}(t)} \left| \sum_{l=1}^{L} \frac{1}{Ld_{\delta}(t)|\mathcal{N}_{ij,\delta}(t)|} y_{ij}^{(l)}(t') - y_{ij}^{*}(t') \right|,$$

clearly $||\bar{\Delta}||_2 \le \max_{i \in [n]} \sum_{j=1}^n |\bar{\Delta}_{ij}(t)| \le \max_{i \in [n]} R_i(t)$ (since $\bar{\Delta}_{ii}(t) = 0$). In order to bound $R_i(t)$ for any given $i \in [n]$, let us first denote

$$S_i = \{ (\xi_{jt'}^{(i)})_{\{i,j\} \in \mathcal{E}_{t'}, t' \in \mathcal{N}_{ij,\delta}(t)} | \xi_{jt'}^{(i)} \in \{-1, 1\} \}.$$
Since the random variables $\frac{\xi_{jt'}^{(i)}}{|\mathcal{N}_{ij,\delta}(t)|} (y_{ij}^{(l)}(t') - y_{ij}^{*}(t'))$ lie almost surely in the interval $[-\frac{1}{N_{\min,\delta}(t)}, \frac{1}{N_{\min,\delta}(t)}]$, we obtain via a simple union bound argument, along with Hoeffding's inequality, that

$$ 
\begin{aligned}
\mathbb{P}(R_i(t) > s) & \le \sum_{\xi^{(i)} \in S_i} \mathbb{P} \left( \sum_{j \ne i} \sum_{t' \in \mathcal{N}_{ij,\delta}(t)} \xi_{jt'}^{(i)} \sum_{l=1}^{L} \frac{1}{Ld_{\delta}(t)|\mathcal{N}_{ij,\delta}(t)|} (y_{ij}^{(l)}(t') - y_{ij}^{*}(t')) > s \right) \\
& \le \sum_{\xi^{(i)} \in S_i} \mathbb{P} \left( \sum_{j \ne i} \sum_{t' \in \mathcal{N}_{ij,\delta}(t)} \xi_{jt'}^{(i)} \sum_{l=1}^{L} \frac{1}{|\mathcal{N}_{ij,\delta}(t)|} (y_{ij}^{(l)}(t') - y_{ij}^{*}(t')) > Ld_{\delta}(t)s \right) \\
& \le \sum_{\xi^{(i)} \in S_i} \exp \left( - \frac{L^2 d_{\delta}^2(t) N_{\min,\delta}^2(t)s^2}{2Ld_{\delta,\max}(t)N_{\max,\delta}(t)} \right) \\
& \le \sum_{\xi^{(i)} \in S_i} \exp \left( - \frac{Ld_{\delta}^2(t) N_{\min,\delta}^2(t)s^2}{2d_{\delta,\max}(t)N_{\max,\delta}(t)} \right).
\end{aligned}
$$

Since $|S_i| \le 2^{d_{\delta,\max}(t)N_{\max,\delta}(t)}$ for each $i \in [n]$, we obtain

$$ \mathbb{P}(R_i(t) > s) \le 2^{d_{\delta,\max}(t)N_{\max,\delta}(t)} \exp \left( - \frac{Ld_{\delta}^2(t) N_{\min,\delta}^2(t)s^2}{2d_{\delta,\max}(t)N_{\max,\delta}(t)} \right) = \exp \left( d_{\delta}(t)N_{\max,\delta}(t) \ln 2 - \frac{Ld_{\delta}^2(t) N_{\min,\delta}^2(t)s^2}{2d_{\delta,\max}(t)N_{\max,\delta}(t)} \right). $$

Applying a union bound over $[n]$ now leads to

$$ \mathbb{P}(\|\bar{\Delta}(t)\|_2 > s) \le n \exp \left( d_{\delta,\max}(t) N_{\max,\delta}(t) \ln 2 - \frac{L d_{\delta}^2(t) N_{\min,\delta}^2(t) s^2}{2 d_{\delta,\max}(t) N_{\max,\delta}(t)} \right). $$

Finally, since $d_{\delta,\max}(t)N_{\max,\delta}(t) \le \log n$ by assumption we obtain for any constant $c_2 > 2$ that

$$ \mathbb{P} \left( \| \bar{\Delta}(t) \|_2 > c_2 \sqrt{ \frac{N_{\max,\delta}(t) d_{\delta,\max}(t) \log n}{L d_{\delta}^2(t) N_{\min,\delta}^2(t)} } \right) \le n^{2-\frac{c_2^2}{2}}. $$

### C.2.3 Bound on $\|\bar{\Delta}(t)\|_2$ when $d_{\delta,\max}(t)N_{\max,\delta}(t) \ge \log n$

In this case, we will use Bernstein's inequality for random matrices (see Theorem 7 in Appendix A). For each $i < j \in [n]$, let $Z_{ij,l}^{(i)}$ be a $n \times n$ matrix with all entries equal to zero except for $Z_{ij,l}^{(i)}$ and $Z_{ji,l}^{(i)}$ defined as follows.

$$ Z_{ij}^{(i,l)}(t') = \frac{C_{ij}^{(i)}(t')}{Ld_{\delta}(t)|\mathcal{N}_{ij,\delta}(t)|} \quad \text{and} \quad Z_{ji}^{(i,l)}(t') = -\frac{C_{ij}^{(i)}(t')}{Ld_{\delta}(t)|\mathcal{N}_{ij,\delta}(t)|} $$

where $C_{ij}^{(i)}(t') = y_{ij}^{(i)}(t') - y_{ij}^{*}(t')$ if $\{i,j\} \in \mathcal{E}_{\delta}(t)$, which is a centered Bernoulli variable, and $C_{ij}^{(i)}(t') = 0$ otherwise. Hence, $\mathbb{E}[Z_{ij,l}^{(i)}] = 0$. Now we can write

$$ \bar{\Delta}(t) = \sum_{\substack{i<j \\ \{i,j\} \in \mathcal{E}_{\delta}(t)}} \sum_{t' \in \mathcal{N}_{ij,\delta}(t)} \sum_{l=1}^{L} Z_{ij,l}^{(i)}(t'). $$
To apply Bernstein's inequality, we first note that $||Z^{ij,l}(t')||_2$ is uniformly bounded for all $t', i, j$.

$$ ||Z^{ij,l}(t')||_2 \le ||Z^{ij,l}(t')||_F = \frac{\sqrt{2}|C_{ij}^{(l)}(t')|}{Ld_{\delta}(t)|\mathcal{N}_{ij,\delta}(t)|} \le \frac{\sqrt{2}}{Ld_{\delta}(t)N_{\min,\delta}(t)} =: B. $$

Since $Z^{ij,l}(t')$ are skew-symmetric and independent matrices, we also have that

$$ \nu = \| \mathbb{E}[\bar{\Delta}(t) \bar{\Delta}(t)^{\top}] \|_2 = \left\| \sum_{i<j} \sum_{t' \in \mathcal{N}_{ij,\delta}(t)} \sum_{l=1}^{L} -\mathbb{E}[Z^{ij,l}(t')^2] \right\|_2. $$

The matrix $(Z^{ij,l}(t'))^2$ only has two non-zero entries, which are at the locations $(i, i)$ and $(j, j)$, and they are equal to $-\frac{C_{ij}^{(l)}(t')^2}{L^2 d_{\delta}^2(t) |\mathcal{N}_{ij,\delta}(t)|^2}$. Then its expectation (which also has only two non-zero entries) satisfies

$$ -\mathbb{E}[Z^{ij,l}(t')^2]_{i,i} = -\mathbb{E}[Z^{ij,l}(t')^2]_{j,j} = \frac{\mathbb{E}[C_{ij}^{(l)}(t')^2]}{L^2 d_{\delta}^2(t) |\mathcal{N}_{ij,\delta}(t)|^2} \le \frac{\mathrm{Var}(C_{ij}^{(l)}(t'))}{L^2 d_{\delta}^2(t) N_{\min,\delta}^2(t)} \le \frac{1}{4L^2 d_{\delta}^2(t) N_{\min,\delta}^2(t)}. $$

since $\mathrm{Var}(C_{ij}^{(l)}(t')) = y_{ij}^*(t')(1 - y_{ij}^*(t')) \le 1/4$. Furthermore, $\mathbb{E}[\bar{\Delta}(t)\bar{\Delta}(t)^{\top}]$ is a diagonal matrix with positive entries, hence

$$ \nu = \max_{k \in [n]} \left| \sum_{i<j} \sum_{t' \in \mathcal{N}_{ij,\delta}(t)} \sum_{l=1}^{L} -\mathbb{E}[Z^{ij}(t')^2] \right|_{kk} \le \frac{Ld_{\delta,\max}(t)N_{\max,\delta}(t)}{4L^2 d_{\delta}^2(t) N_{\min,\delta}^2(t)} = \frac{N_{\max,\delta}(t)d_{\delta,\max}(t)}{4Ld_{\delta}^2(t) N_{\min,\delta}^2(t)}. $$

Finally applying Bernstein's inequality, it holds for any $s > 0$ that

$$ \begin{aligned} \mathbb{P}(\|\bar{\Delta}(t)\|_2 > s) & \le 2n \exp\left(-\frac{3s^2}{6\nu + 2Bs}\right) \\ & \le 2n \exp\left(-\frac{3s^2}{\frac{6N_{\max,\delta}(t)d_{\delta,\max}(t)}{4Ld_{\delta}^2(t)N_{\min,\delta}^2(t)} + \frac{2s\sqrt{2}}{Ld_{\delta}(t)N_{\min,\delta}(t)}}\right) \\ & \le 2n \exp\left(-\frac{6Ld_{\delta}^2(t)N_{\min,\delta}^2(t)s^2}{3N_{\max,\delta}(t)d_{\delta,\max}(t) + 4\sqrt{2}N_{\min,\delta}(t)d_{\delta}(t)s}\right). \end{aligned} $$

Choosing $s = c_3 \sqrt{\frac{N_{\max,\delta}(t)d_{\delta,\max}(t) \log n}{Ld_{\delta}^2(t) N_{\min,\delta}^2(t)}}$ for a constant $c_3 \ge 2$ and since $d_{\delta,\max}(t)N_{\max,\delta}(t) \ge \log n$ by assumption, it holds that

$$ \mathbb{P} \left( \| \bar{\Delta}(t) \|_2 > c_3 \sqrt{\frac{N_{\max,\delta}(t) d_{\delta,\max}(t) \log n}{L d_{\delta}^2(t) N_{\min,\delta}^2(t)}} \right) \le 2n^{1-\frac{6c_3^2}{3+4\sqrt{2}c_3}}. \quad (C.1) $$

Finally, we observe that up to a multiplicative constant, the bounds for $||D(t)||_2$ and $||\bar{\Delta}(t)||_2$ are the same. Thus, for constants $c_1, c_2, c_3 \ge 2$, we have

$$ \|\Delta(t)\|_2 \le (c_1 + \max\{c_2, c_3\}) \sqrt{\frac{N_{\max,\delta}(t) d_{\delta,\max}(t) \log n}{L d_{\delta}^2(t) N_{\min,\delta}^2(t)}}. $$
with probability at least $1 - 2n^{1-\frac{c_1^2}{2}} - \max\{n^{2-\frac{c_2^2}{2}}, 2n^{1-\frac{6c_3^2}{3+4\sqrt{2c_3}}}\}$. In simpler words, there exists a constant $\tilde{C}_1 \ge 15$ such that with probability at least $1 - O(n^{-10})$,

$$ \|\Delta(t)\|_2 \le \tilde{C}_1 \sqrt{\frac{N_{\max,\delta}(t)d_{\delta,\max}(t) \log n}{L d_{\delta}^2(t) N_{\min,\delta}^2(t)}}. \quad (C.2) $$

<figure><img src="image_5.png" /></figure>

### C.3 Proof of Lemma 3

The steps below follow the proof of [23, Lemma 4] and are provided for completeness. We will first show $1 - \lambda_{\max}(\bar{P}(t)) \ge \frac{\xi_{\delta}(t)d_{\delta,\min}(t)}{4d_{\delta}(t)b^3(t)}$ using [23, Lemma 6]. To this end, let us denote $Q_{ij}(t) = \frac{1}{d_i(t)}\mathbf{1}_{\{i,j\} \in \mathcal{E}_{\delta}(t)}$ where $d_i(t)$ denotes the degree of the $i$th vertex in $G_{\delta}(t)$ and $\mathbf{1}$ is an indicator variable. Clearly, $Q(t)$ has a unique stationary distribution as it is irreducible (since $G_{\delta}(t)$ is connected), defined as $\mu_i(t) = \frac{d_i(t)}{\sum_{j=1}^n d_j(t)}$. It satisfies the detailed balanced equations of reversibility

$$ \mu_i(t)Q_{ij}(t) = \mu_j(t)Q_{ji}(t). $$

Using [23, Lemma 6], we obtain

$$ \frac{1 - \lambda_{\max}(\bar{P}(t))}{1 - \lambda_{\max}(Q(t))} \ge \frac{\alpha}{\beta} \quad \text{with } \alpha = \min_{\{i,j\} \in \mathcal{E}_{\delta}(t)} \frac{\pi_i^*(t)\bar{P}_{ij}(t)}{\mu_i(t)Q_{ij}(t)}, \quad \beta = \max_i \frac{\pi_i^*(t)}{\mu_i(t)}. \quad (C.3) $$

Note that the transition matrix $Q$ is equal to the Laplacian $\mathcal{L}_{\delta}(t)$ of the graph $G_{\delta}(t)$ and so $1 - \lambda_{\max}(Q(t)) = \xi_{\delta}(t)(t)$. We now need to suitably bound $\alpha$ and $\beta$.

Denoting $d_i(t)$ to be the degree of node $i$ in $G_{\delta}(t)$, clearly $\sum_{j=1}^n d_j(t) = 2|\mathcal{E}_{\delta}(t)|$. Hence we have

$$ \mu_i(t)Q_{ij}(t) = \frac{1}{\sum_{k=1}^n d_k(t)} \le \frac{1}{|\mathcal{E}_{\delta}(t)|} \quad \text{for } \{i,j\} \in \mathcal{E}_{\delta}(t) \quad \text{and} \quad \mu_i(t) = \frac{d_i(t)}{2|\mathcal{E}_{\delta}(t)|} \ge \frac{d_{\delta,\min}(t)}{2|\mathcal{E}_{\delta}(t)|} \quad \text{for } i \in [n] $$

It is also easy to see that

$$ \pi_i^*(t)\bar{P}_{ij}(t) \ge \frac{1}{2d_{\delta}(t)nb^2(t)} \quad \text{for } \{i,j\} \in \mathcal{E}_{\delta}(t), \quad \text{and} \quad \pi_i^*(t) \le \frac{b(t)}{n} \quad \text{for } i \in [n] $$

thus leading to the bounds

$$ \alpha \ge \frac{|\mathcal{E}_{\delta}(t)|}{2d_{\delta}(t)nb^2(t)} \quad \text{and} \quad \beta \le \frac{2b(t)|\mathcal{E}_{\delta}(t)|}{nd_{\delta,\min}(t)}. $$

Plugging this in (C.3) leads to

$$ 1 - \lambda_{\max}(\bar{P}(t)) \ge \frac{\xi_{\delta}(t)|\mathcal{E}_{\delta}(t)|nd_{\delta,\min}(t)}{4d_{\delta}(t)nb^3(t)|\mathcal{E}_{\delta}(t)|} \ge \frac{\xi_{\delta}(t)d_{\delta,\min}(t)}{4d_{\delta}(t)b^3(t)} > 0 $$

which together with (4.6) readily leads to the stated lower bound on $1 - \rho(t)$.

<figure><img src="image_6.png" /></figure>
### C.4 Proof of Lemma 5

We will bound $\|\pi^*(t)^\top \Delta(t)\|_2$, using the same ideas as in the proof of [7, Theorem 9]. Recall that the entries of $\Delta(t)$ are given by

$$ \Delta_{ij}(t) = \begin{cases} \frac{1}{d_{\delta}(t)|\mathcal{N}_{ij,\delta}(t)|} \sum_{t' \in \mathcal{N}_{ij,\delta}(t)} (y_{ij}(t') - y_{ij}^*(t')) & \text{if } \{i,j\} \in \mathcal{E}_{\delta}(t), \\ -\sum_{k \neq i} \Delta_{ik}(t) & \text{if } i = j, \\ 0 & \text{otherwise.} \end{cases} $$

Let $\Delta^l(t)$ denote the lower triangular part of $\Delta(t)$ and $\Delta^u(t)$ its upper triangular part, both having zeros on the diagonal. Let us define the diagonal matrices

$$ \Delta^{diag,l}(t) = - \mathrm{diag} \left( \sum_{j<i} \Delta_{ij}(t) \right)_{i=1\dots n} \quad \text{and} \quad \Delta^{diag,u}(t) = - \mathrm{diag} \left( \sum_{j>i} \Delta_{ij}(t) \right)_{i=1\dots n} . \quad (C.4) $$

Note that both these matrices have independent entries. The matrix $\Delta(t)$ is now decomposed as

$$ \Delta(t) = \Delta^l(t) + \Delta^u(t) + \Delta^{diag,l}(t) + \Delta^{diag,u}(t) $$

and the triangle inequality implies

$$ \|\pi^*(t)^\top \Delta(t)\|_2 \le \underbrace{\|\pi^*(t)^\top \Delta^l(t)\|_2}_{=I_l} + \underbrace{\|\pi^*(t)^\top \Delta^u(t)\|_2}_{=I_u} + \underbrace{\|\pi^*(t)^\top \Delta^{diag,l}(t)\|_2}_{=I_{diag,l}} + \underbrace{\|\pi^*(t)^\top \Delta^{diag,u}(t)\|_2}_{=I_{diag,u}} \quad (C.5) $$

The goal is to bound each of these terms separately, using the same method for all of them. Let us describe the steps for the first term $I_l = \|\pi^*(t)^\top \Delta^l(t)\|_2$. For all $j \in [n]$,

$$ \left[ \pi^*(t)^\top \Delta^l(t) \right]_j = \sum_{i=1}^{n} \pi_i^*(t) \Delta_{ij}^l(t) = \sum_{i>j} \pi_i^*(t) \Delta_{ij}(t). $$

For any pair $\{i,j\} \in \mathcal{E}_{\delta}(t)$, $\Delta_{ij}(t)$ is a sum of at most $LN_{\max,\delta}(t)$ independent centered random variables, thus $[\pi^*(t)^\top \Delta^l(t)]_j$ is a sum of at most $LN_{\max,\delta}(t) d_j^l(t)$ independent centered random variables, where $d_j^l(t)$ is defined as

$$ d_j^l(t) := |\{i \in [n] \mid \{i,j\} \in \mathcal{E}_{\delta}(t) \text{ and } i > j\}|. $$

By applying Hoeffding's inequality we obtain for any $s > 0$,

$$ \mathbb{P}\left(\left|\left[\pi^*(t)^\top \Delta^l(t)\right]_j\right| > s\right) \le 2 \exp\left(-\frac{2L^2 N_{\min,\delta}^2(t) d_{\delta}^2(t) s^2}{4L N_{\max,\delta}(t) d_j^l(t) \|\pi^*(t)\|_{\infty}^2}\right) \\ \le 2 \exp\left(-\frac{LN_{\min,\delta}^2(t) d_{\delta}^2(t) s^2}{2N_{\max,\delta}(t) d_j^l(t) \|\pi^*(t)\|_{\infty}^2}\right). $$

Hence, $[\pi^*(t)^\top \Delta^l(t)]_j$ can be seen as a sub-Gaussian random variable with sub-Gaussian norm at most $c_1 \sigma_j$ with $c_1 > 0$ a constant, and $\sigma_j^2 = \frac{2d_j^l(t) N_{\max,\delta}(t) \|\pi^*(t)\|_{\infty}^2}{L d_{\delta}^2(t) N_{\min,\delta}^2(t)}$. As $d_j^l(t) < d_{\delta,\max}(t)$, it holds that

$$ \sigma_j^2 \le 2 \frac{N_{\max,\delta}(t) d_{\delta,\max}(t)}{L d_{\delta}^2(t) N_{\min,\delta}^2(t)} \|\pi^*(t)\|_{\infty}^2 $$
Denoting $\sigma^2 = 2c_1 \frac{N_{\max,\delta}(t)d_{\delta,\max}(t)}{Ld_{\delta}^2(t)N_{\min,\delta}^2(t)} \|\pi^*(t)\|_{\infty}^2$, then we see that $\pi^*(t)^{\top}\Delta^l(t)$ is a random vector with independent, centered $\sigma$- sub-gaussian entries. Then,

$$I_l^2 = \|\pi^*(t)^{\top}\Delta^l(t)\|_2^2 = \sum_{j=1}^{n} [\pi^*(t)^{\top}\Delta^l(t)]_j^2$$

is a quadratic form of a sub-gaussian vector, so its expectation satisfies $\mathbb{E}[I_l^2] \le c_2 n\sigma^2$ for some constant $c_2 > 0$.

Using the Hanson-Wright inequality [26, Theorem 1.1], we have for any $s > 0$ that

$$\mathbb{P}(|I_l^2 - \mathbb{E}[I_l^2]| > s) \le 2 \exp \left( -c_3 \min \left( \frac{s^2}{\sigma^4 n}, \frac{s}{\sigma^2} \right) \right) = 2 \exp \left( -c_3 \frac{s}{\sigma^2} \min \left( \frac{s}{\sigma^2 n}, 1 \right) \right)$$

for $c_3 > 0$ a constant. Choosing $s = \sqrt{\frac{C_1}{c_3}}\sigma^2\sqrt{n\log n}$ with $C_1 \ge 1$ a constant, note that $\frac{s}{n\sigma^2} = \sqrt{\frac{C_1}{c_3}}\sqrt{\frac{\log n}{n}} \le 1$ if $n \ge \frac{C_1}{c_3}\log n$, and so,

$$\mathbb{P}\left( |I_l^2 - \mathbb{E}[I_l^2]| > \sqrt{\frac{C_1}{c_3}}\sigma^2\sqrt{n\log n} \right) \le 2\exp^{-C_1\log n} = 2n^{-C_1}.$$

Then, it holds with probability at least $1 - O(n^{-C_1})$ that

$$\begin{aligned}
I_l^2 &\le \mathbb{E}[I_l^2] + \sigma^2 \sqrt{n \log n} \\
&\le c_2 n \sigma^2 + \sqrt{\frac{C_1}{c_3}} \sigma^2 \sqrt{n \log n} \\
&\le n \sigma^2 (c_2 + \sqrt{\frac{C_1}{c_3}} \sqrt{\frac{\log n}{n}}) \\
&\le C_2 \frac{n N_{\max,\delta}(t) d_{\delta,\max}(t)}{L d_{\delta}^2(t) N_{\min,\delta}^2(t)} \|\pi^*(t)\|_{\infty}^2 \\
&\le C_2 \frac{N_{\max,\delta}(t) d_{\delta,\max}(t) b^2(t)}{L d_{\delta}^2(t) N_{\min,\delta}^2(t)} \|\pi^*(t)\|_2^2
\end{aligned}$$

$$\|\pi^*(t)\|_2^2 \ge n \pi_{\min}^*(t)^2 \ge \frac{n}{b^2(t)} \pi_{\max}^*(t)^2 = \frac{n}{b^2(t)} \|\pi^*(t)\|_{\infty}^2.$$

The same analysis leads to the same bound for $I_u^2, I_{\text{diag},l}^2, I_{\text{diag},u}^2$. Thus for any constant $C_1 \ge 1$, there exists a constant $C_2 \ge 1$ such that it holds with probability greater than $1 - O(n^{-C_1})$ that

$$\|\pi^*(t)^{\top}\Delta(t)\|_2 \le C_2 b(t) \sqrt{\frac{N_{\max,\delta}(t)d_{\delta,\max}(t)}{Ld_{\delta}^2(t)N_{\min,\delta}^2(t)}} \|\pi^*(t)\|_2. \quad (C.6)$$

# D Proofs of results in Section 5

The proofs are largely inspired from the proofs in [7, Section C], which are adapted to our problem setup.
## D.1 Proof of Lemma 6

Consider any given collection of graphs $(G_{t'})_{t' \in \mathcal{N}_\delta(t)}$. Then following the same ideas as in [7, Lemma 2], we have for any given $m \in [n]$ that

$$ 
\begin{aligned}
I_0^m &= \sum_{j=1}^n \pi_j^*(t)(\hat{P}_{jm}(t) - \hat{P}_{jm}^*(t)) \\
&= \sum_{j:j \neq m} \pi_j^*(t)(\hat{P}_{jm}(t) - \hat{P}_{jm}^*(t)) + \pi_m^*(t)(\hat{P}_{mm}(t) - \hat{P}_{mm}^*(t)) \\
&= \sum_{j:j \neq m} (\pi_j^*(t) + \pi_m^*(t))(\hat{P}_{jm}(t) - \hat{P}_{jm}^*(t)) \\
&= \sum_{j:j \neq m} (\pi_j^*(t) + \pi_m^*(t)) \sum_{t' \in \mathcal{N}_{jm,\delta}(t)} \sum_{l=1}^L \left( \frac{y_{jm}^{(l)}(t') - y_{jm}^*(t')}{Ld_\delta(t)|\mathcal{N}_{jm,\delta}(t)|} \right).
\end{aligned}
$$

Applying Hoeffding's inequality, it holds for any constant $C_1 \ge 2$,

$$ \mathbb{P}(|I_0^m| > C_1 \sqrt{\frac{N_{\max,\delta}(t)d_{\delta,\max}(t) \log n}{Ld_\delta^2(t)N_{\min,\delta}^2(t)}} \|\pi^*(t)\|_\infty) \le 2n^{-\frac{C_1^2}{2}}. $$

If $G_{t'} \sim \mathcal{G}(n, p(t'))$ for all $t' \in \mathcal{N}_\delta(t)$, then using Lemma 11, it holds with probability at least $1 - O(n^{-10})$ that $\frac{N_{\max,\delta}(t)d_{\delta,\max}(t)}{d_\delta^2(t)N_{\min,\delta}^2(t)} \le \frac{4}{3np_\delta(t)p_{\delta,\text{sum}}(t)}$. Hence after a union bound over $[n]$, we conclude that there exists a suitably large constant $C_2 \ge 1$ such that

$$ \mathbb{P} \left( \forall m \in [n] : |I_0^m| > C_2 \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} \|\pi^*(t)\|_\infty \right) \le O(n^{-9}). $$

□

## D.2 Proof of Lemma 7

For any given collection of graphs $(G_{t'})_{t' \in \mathcal{N}_\delta(t)}$, recall that the entries of $\hat{P}^*(t) - \bar{P}(t)$ have already been bounded in the proof of Lemma 1 as follows

$$ \left| \hat{P}_{ij}^*(t) - \bar{P}_{ij}(t) \right| \le \begin{cases} \frac{M\delta}{Td_\delta(t)} & \text{if } i \neq j, \{i,j\} \in \mathcal{E}_\delta(t) \\ \frac{M\delta}{T} & \text{if } i = j \\ 0 & \text{otherwise.} \end{cases} \quad (D.1) $$

This implies for all $m \in [n]$ that

$$ |I_1^m| = \left| \sum_{j:j \neq m} \pi_j^*(t)(\hat{P}_{jm}^*(t) - \bar{P}_{jm}(t)) \right| \le \|\pi^*(t)\|_\infty \left( d_{\delta,\max}(t) \frac{M\delta}{Td_\delta(t)} + \frac{M\delta}{T} \right). $$

Due to Lemma 11, we know that $3np_\delta(t) = d_\delta(t) \ge d_{\delta,\max}(t)$ holds with probability at least $1 - O(n^{-10})$. This leads to the statement in Lemma 7.

□
### D.3 Proof of Lemma 8

We will follow the steps in [7, Lemma 3]. Recalling that $I_2^m = (\hat{\pi}_m(t) - \pi_m^*(t)) \hat{P}_{mm}(t)$, let us bound $\hat{P}_{mm}(t)$ using the decomposition

$$ \hat{P}_{mm}(t) = (\hat{P}_{mm}(t) - \hat{P}_{mm}^*(t)) + \hat{P}_{mm}^*(t). $$

For any given collection of graphs $(G_{t'})_{t' \in \mathcal{N}_\delta(t)}$, note that

$$ \hat{P}_{mm}(t) - \hat{P}_{mm}^*(t) = -\frac{1}{Ld_\delta(t)} \sum_{j: j \neq m} \sum_{t' \in \mathcal{N}_{jm, \delta}(t)} \sum_{l=1}^{L} \frac{1}{|\mathcal{N}_{jm, \delta}(t)|} (y_{mj}^{(l)}(t') - y_{mj}^*(t')) $$

is a sum of independent, centered random variables. Then Hoeffding's inequality implies for each $m \in [n]$,

$$ \mathbb{P}\left(\left|\hat{P}_{mm}(t) - \hat{P}_{mm}^*(t)\right| > C_1 \sqrt{\frac{N_{\max,\delta}(t)d_{\delta,\max}(t) \log n}{Ld_\delta^2(t)N_{\min,\delta}^2(t)}}\right) \le 2n^{-\frac{C_1^2}{2}} $$

for any constant $C_1 \ge 2$.

If $G_{t'} \sim \mathcal{G}(n, p(t'))$ for all $t' \in \mathcal{N}_\delta(t)$, then once again combining this result with Lemma 11, we see (after a union bound over $[n]$) that there exists a suitably large constant $C_2 \ge 1$ such that with probability at least $1 - O(n^{-9})$,

$$ \left| \hat{P}_{mm}(t) - \hat{P}_{mm}^*(t) \right| \le C_2 \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}}, \quad \forall m \in [n]. $$

Now we bound $|\hat{P}_{mm}^*(t)|$. Recalling that $b_{\max,\delta}(t) = \max_{t' \in \mathcal{N}_\delta(t)} b(t')$ and $b(t') = \max_{i,j \in [n]} \frac{w_{t',i}^*}{w_{t',j}^*}$, we have for any given collection of graphs $(G_{t'})_{t' \in \mathcal{N}_\delta(t)}$ that all $m \in [n]$,

$$ \begin{aligned} \hat{P}_{mm}^*(t) &= 1 - \sum_{j: j \neq m} \frac{1}{d_\delta(t)|\mathcal{N}_{jm, \delta}(t)|} \sum_{t' \in \mathcal{N}_{jm, \delta}(t)} y_{mj}^*(t') \\ &= 1 - \sum_{j: j \neq m} \sum_{t' \in \mathcal{N}_{jm, \delta}(t)} \frac{1}{d_\delta(t)|\mathcal{N}_{jm, \delta}(t)|} \left( \frac{w_{t',j}^*}{w_{t',j}^* + w_{t',m}^*} \right) \\ &\le 1 - \sum_{j: j \neq m} \sum_{t' \in \mathcal{N}_{jm, \delta}(t)} \frac{1}{d_\delta(t)|\mathcal{N}_{jm, \delta}(t)|} \left( \frac{1}{1 + b(t')} \right) \\ &\le 1 - \frac{d_{\delta, \min}(t)}{d_\delta(t)} \frac{1}{2b_{\max, \delta}(t)}. \end{aligned} $$

If $G_{t'} \sim \mathcal{G}(n, p(t'))$ for all $t' \in \mathcal{N}_\delta(t)$, then we have using Lemma 11 that $d_{\delta, \min}(t) \ge np_\delta(t)/2$ with probability at least $1 - O(n^{-10})$, this in turn implies that $|\hat{P}_{mm}^*(t)| \le 1 - \frac{1}{12b_{\max, \delta}(t)}$ for all $m \in [n]$.

Hence, we conclude that there exists a suitably large constant $C_2 \ge 1$ such that with probability at least $1 - O(n^{-9})$,

$$ |I_2^m| \le \left( 1 - \frac{1}{12b_{\max,\delta}(t)} + C_2 \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} \right) \|\hat{\pi}(t) - \pi^*(t)\|_{\infty}, \quad \forall m \in [n]. $$

□
## D.4 Proof of Lemma 9

This proof follows the lines of the proof of [7, Lemma 4]. Consider any given collection of graphs $(G_{t'})_{t' \in \mathcal{N}_{\delta}(t)}$. Using the properties of $\|\cdot\|_{\pi^{*}(t)}$ in (2.1) and Theorem 4 with the matrices $\hat{P}(t), \hat{P}^{(m)}(t), \bar{P}(t)$, it holds for all $m \in [n]$ that

$$ \begin{aligned} \|\hat{\pi}(t) - \hat{\pi}^{(m)}(t)\|_{2} &\le \frac{1}{\sqrt{\pi_{\min}^{*}(t)}} \left( \frac{1}{1 - \lambda_{\max}(\bar{P}(t)) - \|\hat{P}(t) - \bar{P}(t)\|_{\pi^{*}(t)}} \right) \|\hat{\pi}^{(m)}(t)^{\top}(\hat{P}^{(m)}(t) - \hat{P}(t))\|_{\pi^{*}(t)} \\ &\le \sqrt{\frac{\pi_{\max}^{*}(t)}{\pi_{\min}^{*}(t)}} \left( \frac{8d_{\delta}(t)b^{3}(t)}{\xi_{\delta}(t)d_{\delta,\min}(t)} \right) \|\hat{\pi}^{(m)}(t)^{\top}(\hat{P}^{(m)}(t) - \hat{P}(t))\|_{2} \end{aligned} \quad (D.2) $$

where we used (4.7).

In order to bound $\|\hat{\pi}^{(m)}(t)^{\top}(\hat{P}^{(m)}(t) - \hat{P}(t))\|_{2}$, we need to get rid of the statistical dependency between $\hat{P}^{(m)}(t)$ and $\hat{\pi}^{(m)}(t)$. To do so, let us introduce $\hat{P}^{(m),\mathcal{G}}(t)$, defined for $i \ne j$ as

$$ \hat{P}_{ij}^{(m),\mathcal{G}}(t) = \begin{cases} \hat{P}_{ij}(t) & \text{if } i \ne m, j \ne m \\ \frac{1}{d_{\delta}(t)} \frac{w_{t,j}^{*} w_{t,i}^{*}}{w_{t,j}^{*} + w_{t,i}^{*}} \mathbf{1}_{\{i,j\} \in \mathcal{E}_{\delta}(t)} & \text{if } i = m \text{ or } j = m. \end{cases} $$

Its diagonal entries are defined as

$$ \hat{P}_{ii}^{(m),\mathcal{G}}(t) = 1 - \sum_{k:k \ne i} \hat{P}_{ik}^{(m),\mathcal{G}}(t) $$

to ensure that $\hat{P}^{(m),\mathcal{G}}(t)$ is a transition matrix. Note that the $m^{th}$ line and columns of $\hat{P}^{(m)}(t)$ are the expectations of those in $\hat{P}^{(m),\mathcal{G}}(t)$. Moreover, $\hat{P}_{\cdot m}^{(m),\mathcal{G}}(t) = \bar{P}_{\cdot m}(t)$ and $\hat{P}_{\cdot m}^{(m),\mathcal{G}}(t) = \bar{P}_{\cdot m}(t)$. Now starting with the decomposition

$$ \hat{\pi}^{(m)}(t)^{\top}(\hat{P}^{(m)}(t) - \hat{P}(t)) = \underbrace{\hat{\pi}^{(m)}(t)^{\top}(\hat{P}(t) - \hat{P}^{(m),\mathcal{G}}(t))}_{J_1^m} + \underbrace{\hat{\pi}^{(m)}(t)^{\top}(\hat{P}^{(m),\mathcal{G}}(t) - \hat{P}^{(m)}(t))}_{J_2^m} $$

we will separately bound $\|J_1^m\|_2, \|J_2^m\|_2$.

### D.4.1 Bound on $\|J_1^m\|_2$

Consider any given collection of graphs $(G_{t'})_{t' \in \mathcal{N}_{\delta}(t)}$. Since $\hat{P}_{\cdot m}^{(m),\mathcal{G}}(t) = \bar{P}_{\cdot m}(t)$, we have the decomposition

$$ \begin{aligned} J_{1,m}^{m} &= \sum_{j=1}^{n} \hat{\pi}_{j}^{(m)}(t) (\hat{P}_{jm}(t) - \bar{P}_{jm}(t)) \\ &= \sum_{j:j \ne m} \hat{\pi}_{j}^{(m)}(t) (\hat{P}_{jm}(t) - \bar{P}_{jm}(t)) + \hat{\pi}_{m}^{(m)}(t) (\hat{P}_{mm}(t) - \bar{P}_{mm}(t)) \\ &= \sum_{j:j \ne m} (\hat{\pi}_{j}^{(m)}(t) + \hat{\pi}_{m}^{(m)}(t)) (\hat{P}_{jm}(t) - \bar{P}_{jm}(t)) \\ &= \sum_{j:j \ne m} (\hat{\pi}_{j}^{(m)}(t) + \hat{\pi}_{m}^{(m)}(t)) (\hat{P}_{jm}(t) - \hat{P}_{jm}^{*}(t)) + \sum_{j:j \ne m} (\hat{\pi}_{j}^{(m)}(t) + \hat{\pi}_{m}^{(m)}(t)) (\hat{P}_{jm}^{*}(t) - \bar{P}_{jm}(t)). \end{aligned} $$
Since $\hat{\pi}^{(m)}(t)$ is independent of $(\hat{P}_{jm}(t) - \hat{P}_{jm}^{*}(t))_{j \neq m}$, the first term can be bounded using Hoeffding's inequality as

$$ \mathbb{P}\left(\left|\sum_{j:j \neq m} (\hat{\pi}_j^{(m)}(t) + \hat{\pi}_m^{(m)}(t)) (\hat{P}_{jm}(t) - \hat{P}_{jm}^{*}(t))\right| > C_1 \|\hat{\pi}^{(m)}(t)\|_{\infty} \sqrt{\frac{N_{\max, \delta}(t) d_{\delta, \max}(t) \log n}{L d_{\delta}^2(t) N_{\min, \delta}^2(t)}}\right) \le 2n^{-\frac{C_1^2}{2}} $$

for any constant $C_1 \ge 2$. Moreover, due to (D.1), we have for all $m$ that

$$ \sum_{j:j \neq m} (\hat{\pi}_j^{(m)}(t) + \hat{\pi}_m^{(m)}(t)) (\hat{P}_{jm}^{*}(t) - \bar{P}_{jm}(t)) \le 2d_{\delta, \max}(t) \left(\frac{M\delta}{Td_{\delta}(t)}\right) \|\hat{\pi}^{(m)}(t)\|_{\infty}. $$

Hence for any given $m \in [n]$, it holds with probability at least $1 - O(n^{-\frac{C_1^2}{2}})$ that

$$ |J_{1,m}^m| \le \left( \frac{2d_{\delta,\max}(t)M\delta}{Td_{\delta}(t)} + C_1 \sqrt{\frac{N_{\max,\delta}(t)d_{\delta,\max}(t)\log n}{Ld_{\delta}^2(t)N_{\min,\delta}^2(t)}} \right) \|\hat{\pi}^{(m)}(t)\|_{\infty}. \quad (D.3) $$

Let us now bound the other coefficients of $J_1^m$. For any $j \neq m$,

$$ J_{1,j}^m = \sum_{i:i \neq j} \hat{\pi}_i^{(m)}(t) (\hat{P}_{ij}(t) - \hat{P}_{ij}^{(m),\mathcal{G}}(t)) + \hat{\pi}_j^{(m)}(t) (\hat{P}_{jj}(t) - \hat{P}_{jj}^{(m),\mathcal{G}}(t)) $$

$$ = \hat{\pi}_m^{(m)}(t) (\hat{P}_{mj}(t) - \hat{P}_{mj}^{(m),\mathcal{G}}(t)) + \hat{\pi}_j^{(m)}(t) (\hat{P}_{jj}(t) - \hat{P}_{jj}^{(m),\mathcal{G}}(t)) $$

$$ = (\hat{\pi}_m^{(m)}(t) + \hat{\pi}_j^{(m)}(t)) (\hat{P}_{mj}(t) - \hat{P}_{mj}^{*}(t) + \hat{P}_{mj}^{*}(t) - \bar{P}_{mj}(t)) $$

$$ \implies |J_{1,j}^m| \le 2 \|\hat{\pi}^{(m)}(t)\|_{\infty} \left( |\hat{P}_{mj}(t) - \hat{P}_{mj}^{*}(t)| + \frac{M\delta}{Td_{\delta}(t)} \right) $$

where we used (D.1) and the fact that the $m^{th}$ rows of $\hat{P}^{(m),\mathcal{G}}(t)$ and $\bar{P}(t)$ are identical. Using Hoeffding's inequality and a union bound, we have

$$ \mathbb{P}\left(\max_{j \neq m} |\hat{P}_{mj}(t) - \hat{P}_{mj}^{*}(t)| > C_2 \sqrt{\frac{N_{\max, \delta}(t) \log n}{L d_{\delta}^2(t) N_{\min, \delta}^2(t)}}\right) \le 2n^{1-\frac{C_2^2}{2}} $$

for any constant $C_2 \ge 2$. Hence, for any $j \neq m$, we have with probability at least $1 - O(n^{1-\frac{C_2^2}{2}})$ that

$$ |J_{1,j}^m| \le \begin{cases} 2 \|\hat{\pi}^{(m)}(t)\|_{\infty} \left( \frac{M\delta}{Td_{\delta}(t)} + C_2 \sqrt{\frac{N_{\max,\delta}(t) \log n}{Ld_{\delta}^2(t) N_{\min,\delta}^2(t)}} \right) & \text{if } \{j, m\} \in \mathcal{E}_{\delta}(t) \\ 0 & \text{otherwise.} \end{cases} \quad (D.4) $$

Finally combining (D.3), (D.4) along with Lemma 11, we see after taking a union bound over $[n]$ that when $G_{t'} \sim \mathcal{G}(n, p(t'))$ for all $t' \in \mathcal{N}_{\delta}(t)$, then there exists a suitably large constant $C_3 \ge 1$ such that with probability at least $1 - O(n^{-9})$,

$$ \|J_1^m\|_2 \le \left( \frac{4M\delta}{T} + C_3 \sqrt{\frac{\log n}{Lnp_{\delta}(t)p_{\delta,\text{sum}}(t)}} \right) \|\hat{\pi}^{(m)}(t)\|_{\infty} \quad \forall m \in [n]. \quad (D.5) $$
### D.4.2 Bound on $||J_2^m||_2$

Consider again any given collection of graphs $(G_{t'})_{t' \in \mathcal{N}_\delta(t)}$. Note that the off-diagonal entries of $\hat{P}^{(m)}(t) - \hat{P}^{(m),\mathcal{G}}(t)$ are non zero if and only if they belong to the $m^{th}$ row /column. Then for $i \neq j$ with $i = m$ or $j = m$, we have

$$ (\hat{P}^{(m)}(t) - \hat{P}^{(m),\mathcal{G}}(t))_{ij} = \frac{1}{d_{\delta}(t)} (p_{\delta}(t) - \mathbf{1}_{\{i,j\} \in \mathcal{E}_{\delta}(t)}) y_{ij}^{*}(t) $$

Moreover, the diagonal entries are given as

$$ (\hat{P}^{(m)}(t) - \hat{P}^{(m),\mathcal{G}}(t))_{ii} = \begin{cases} -\frac{1}{d_{\delta}(t)} \sum_{j:j \neq m} (p_{\delta}(t) - \mathbf{1}_{\{m,j\} \in \mathcal{E}_{\delta}(t)}) y_{mj}^{*}(t) & \text{if } i = m \\ -\frac{1}{d_{\delta}(t)} (p_{\delta}(t) - \mathbf{1}_{\{i,m\} \in \mathcal{E}_{\delta}(t)}) y_{im}^{*}(t) & \text{if } i \neq m. \end{cases} $$

Let us first show that $\pi^*(t)^\top (\hat{P}^{(m)}(t) - \hat{P}^{(m),\mathcal{G}}(t)) = 0$ by considering the cases below.

1. For any $j \neq m$,

$$ \begin{aligned} [\pi^*(t)^\top (\hat{P}^{(m)}(t) - \hat{P}^{(m),\mathcal{G}}(t))]_j &= \frac{\pi_m^*(t)}{d_\delta(t)} (p_\delta(t) - \mathbf{1}_{\{m,j\} \in \mathcal{E}_\delta(t)}) y_{mj}^*(t) - \frac{\pi_j^*(t)}{d_\delta(t)} (p_\delta(t) - \mathbf{1}_{\{j,m\} \in \mathcal{E}_\delta(t)}) y_{jm}^*(t) \\ &= \frac{p_\delta(t) - \mathbf{1}_{\{m,j\} \in \mathcal{E}_\delta(t)}}{d_\delta(t)} (\pi_m^*(t) y_{mj}^*(t) - \pi_j^*(t) y_{jm}^*(t)) \\ &= 0. \end{aligned} $$

This equality comes from the definition of $\pi^*(t)$ and $y_{jm}^*(t)$.

2. For $j = m$,

$$ \begin{aligned} & [\pi^*(t)^\top (\hat{P}^{(m)}(t) - \hat{P}^{(m),\mathcal{G}}(t))]_m \\ &= -\frac{\pi_m^*(t)}{d_\delta(t)} \sum_{k:k \neq m} (p_\delta(t) - \mathbf{1}_{\{m,k\} \in \mathcal{E}_\delta(t)}) y_{mk}^*(t) + \frac{1}{d_\delta(t)} \sum_{k:k \neq m} \pi_k^*(t) (p_\delta(t) - \mathbf{1}_{\{k,m\} \in \mathcal{E}_\delta(t)}) y_{km}^*(t) \\ &= \frac{1}{d_\delta(t)} \sum_{k:k \neq m} (p_\delta(t) - \mathbf{1}_{\{k,m\} \in \mathcal{E}_\delta(t)}) (\pi_k^*(t) y_{mk}^*(t) - \pi_m^*(t) y_{mk}^*(t)) \\ &= 0. \end{aligned} $$

Hence, $\pi^*(t)^\top (\hat{P}^{(m)}(t) - \hat{P}^{(m),\mathcal{G}}(t)) = 0$ and it holds that

$$ J_2^m = (\hat{\pi}^{(m)}(t) - \pi^*(t))^\top (\hat{P}^{(m)}(t) - \hat{P}^{(m),\mathcal{G}}(t)). $$

Now for $j \neq m$,

$$ \begin{aligned} |J_{2,j}^m| &= \left| \left(\hat{\pi}_m^{(m)}(t) - \pi_m^*(t)\right) \left(\frac{p_\delta(t) - \mathbf{1}_{\{m,j\} \in \mathcal{E}_\delta(t)}}{d_\delta(t)}\right) y_{mj}^*(t) - \left(\hat{\pi}_j^{(m)}(t) - \pi_j^*(t)\right) \left(\frac{p_\delta(t) - \mathbf{1}_{\{j,m\} \in \mathcal{E}_\delta(t)}}{d_\delta(t)}\right) y_{jm}^*(t) \right| \\ &\le \begin{cases} \frac{4}{d_\delta(t)} \|\hat{\pi}^{(m)}(t) - \pi^*(t)\|_{\infty} & \text{if } \{j,m\} \in \mathcal{E}_\delta(t) \\ 2 \frac{p_\delta(t)}{d_\delta(t)} \|\hat{\pi}^{(m)}(t) - \pi^*(t)\|_{\infty} & \text{if } \{j,m\} \notin \mathcal{E}_\delta(t). \end{cases} \end{aligned} $$
Moreover, we can bound $|J_{2,m}^m|$ as follows.

$$|J_{2,m}^m| = \left| - \left( \hat{\pi}_m^{(m)}(t) - \pi_m^*(t) \right) \sum_{j:j \neq m} \frac{1}{d_\delta(t)} (p_\delta(t) - \mathbf{1}_{\{m,j\} \in \mathcal{E}_\delta(t)}) y_{mj}^*(t) \right. \\ \left. + \sum_{j:j \neq m} \left( \hat{\pi}_j^{(m)}(t) - \pi_j^*(t) \right) \frac{1}{d_\delta(t)} (p_\delta(t) - \mathbf{1}_{\{j,m\} \in \mathcal{E}_\delta(t)}) y_{jm}^*(t) \right| \\ \le |J_3^m| + |J_4^m|$$

where

$$J_3^m = \sum_{j:j \neq m} \left( \hat{\pi}_m^{(m)}(t) - \pi_m^*(t) \right) \frac{1}{d_\delta(t)} (p_\delta(t) - \mathbf{1}_{\{m,j\} \in \mathcal{E}_\delta(t)}) y_{mj}^*(t),$$

and

$$J_4^m = \sum_{j:j \neq m} \left( \hat{\pi}_j^{(m)}(t) - \pi_j^*(t) \right) \frac{1}{d_\delta(t)} (p_\delta(t) - \mathbf{1}_{\{j,m\} \in \mathcal{E}_\delta(t)}) y_{jm}^*(t).$$

When $G_{t'} \sim \mathcal{G}(n, p(t'))$ for all $t' \in \mathcal{N}_\delta(t)$, then denoting $\xi_j^{(m)} = (\hat{\pi}_m^{(m)}(t) - \pi_m^*(t)) \frac{y_{mj}^*(t)}{d_\delta(t)}$, a bound on $J_3^m$ can be given by Bernstein's inequality, introducing the random variables $Z_j^{(m)} = \xi_j^{(m)}(p_\delta(t) - \mathbf{1}_{\{j,m\} \in \mathcal{E}_\delta(t)})$. Since for all $j \neq m$,

$$|Z_j^{(m)}| \le \frac{1}{d_\delta(t)} \|\hat{\pi}^{(m)}(t) - \pi^*(t)\|_{\infty}, \quad \mathbb{E}[Z_j^{(m)}]^2 \le \frac{\|\hat{\pi}^{(m)}(t) - \pi^*(t)\|_{\infty}^2}{d_\delta^2(t)} p_\delta(t)$$

hence there exists a constant $c \ge 2$ such that it holds with probability at least $1 - 2n^{-3c/2}$ that

$$|J_3^m| \le c \frac{\sqrt{np_\delta(t) \log n} + \log n}{d_\delta(t)} \|\hat{\pi}^{(m)}(t) - \pi^*(t)\|_{\infty}. \quad (D.6)$$

The same bound holds for $|J_4^m|$. Using Lemma 11, we then see that there exists a suitably large constant $C_4 \ge 1$ such that with probability at least $1 - O(n^{-10})$,

$$\|J_2^m\|_2 \le \left( C_4 \frac{\sqrt{np_\delta(t) \log n} + \log n}{d_\delta(t)} + \frac{4\sqrt{d_\delta(t)}}{d_\delta(t)} + \frac{2p_\delta(t)\sqrt{n}}{d_\delta(t)} \right) \|\hat{\pi}^{(m)}(t) - \pi^*(t)\|_{\infty}.$$

Taking the union bound over $[n]$ and performing a minor simplification of the previous bound, we finally observe that with probability at least $1 - O(n^{-9})$,

$$\|J_2^m\|_2 \le C_5 \sqrt{\frac{\log n}{np_\delta(t)}} \|\hat{\pi}^{(m)}(t) - \pi^*(t)\|_{\infty} \quad \forall m \in [n]. \quad (D.7)$$

for some suitably large constant $C_5 \ge 1$.

### D.4.3 Putting it together

Using Lemma 11, we know that with probability at least $1 - O(n^{-10})$,

$$\sqrt{\frac{\pi_{\max}^*(t)}{\pi_{\min}^*(t)}} \left( \frac{8d_\delta(t)b^3(t)}{\xi_\delta(t)d_{\delta,\min}(t)} \right) \le 96b^{\frac{7}{2}}(t). \quad (D.8)$$
Finally, combining (D.2), (D.5), (D.7) and (D.8), we have that with probability at least $1 - O(n^{-9})$, it holds for all $m \in [n]$,

$$
\begin{aligned}
\|\hat{\pi}^{(m)}(t) - \hat{\pi}(t)\|_2 & \le 96b^{\frac{5}{2}}(t)(\|J_1^m\|_2 + \|J_2^m\|_2) \\
& \le 96b^{\frac{5}{2}}(t) \left( \frac{4M\delta}{T} + C_3 \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} \right) \|\hat{\pi}^{(m)}(t)\|_\infty \\
& \phantom{\le 96b^{\frac{5}{2}}(t) \left(} + 96C_5 b^{\frac{5}{2}}(t) \sqrt{\frac{\log n}{np_\delta(t)}} \|\hat{\pi}^{(m)}(t) - \pi^*(t)\|_\infty \\
& \le 96b^{\frac{5}{2}}(t) \left( \frac{4M\delta}{T} + C_3 \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} \right) \|\pi^*(t)\|_\infty \\
& \phantom{\le 96b^{\frac{5}{2}}(t) \left(} + 96b^{\frac{5}{2}}(t) \left( \frac{4M\delta}{T} + C_3 \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} + C_5 \sqrt{\frac{\log n}{np_\delta(t)}} \right) \|\hat{\pi}^{(m)}(t) - \pi^*(t)\|_\infty \\
& \le 96b^{\frac{5}{2}}(t) \left( \frac{4M\delta}{T} + C_3 \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} \right) \|\pi^*(t)\|_\infty \\
& \phantom{\le 96b^{\frac{5}{2}}(t) \left(} + 96b^{\frac{5}{2}}(t) \left( \frac{4M\delta}{T} + C_6 \sqrt{\frac{\log n}{np_\delta(t)}} \right) \|\hat{\pi}^{(m)}(t) - \pi^*(t)\|_\infty
\end{aligned}
$$

for some suitably large constant $C_6 \ge 1$. Condition (5.4) further implies that for all $m \in [n]$,

$$ \|\hat{\pi}^{(m)}(t) - \hat{\pi}(t)\|_2 \le 96b^{\frac{5}{2}}(t) \left( \frac{4M\delta}{T} + C_3 \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} \right) \|\pi^*(t)\|_\infty + \frac{1}{2} \|\hat{\pi}^{(m)}(t) - \pi^*(t)\|_\infty. $$

Finally, the triangular inequality

$$ \|\hat{\pi}^{(m)}(t) - \pi^*(t)\|_{\infty} \le \|\hat{\pi}^{(m)}(t) - \hat{\pi}(t)\|_2 + \|\hat{\pi}(t) - \pi^*(t)\|_{\infty} $$

implies that

$$
\begin{aligned}
\|\hat{\pi}^{(m)}(t) - \hat{\pi}(t)\|_2 & \le 96b^{\frac{5}{2}}(t) \left( \frac{4M\delta}{T} + C_3 \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} \right) \|\pi^*(t)\|_\infty + \frac{1}{2} \left( \|\hat{\pi}^{(m)}(t) - \hat{\pi}(t)\|_2 + \|\hat{\pi}(t) - \pi^*(t)\|_\infty \right) \\
& \le 192b^{\frac{5}{2}}(t) \left( \frac{4M\delta}{T} + C_3 \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} \right) \|\pi^*(t)\|_\infty + \|\hat{\pi}(t) - \pi^*(t)\|_\infty.
\end{aligned}
$$

## D.5 Proof of Lemma 10

Since we follow the same steps as in [6], let us introduce some new notations in our particular setup. For all $t \in \mathcal{T}$, all $i \ne j$ and all $l \in [L]$, let us denote

$$ \tilde{y}_{ij}^{(l)}(t) \sim \mathcal{B}(y_{ij}^{*}(t)) \quad \text{and} \quad \tilde{y}_{ij}(t) = \frac{1}{L} \sum_{l=1}^{L} \tilde{y}_{ij}^{(l)}(t). \qquad (\text{D.9}) $$

so for all $t' \in \mathcal{N}_\delta(t)$, we have

$$ y_{ij}(t') = \tilde{y}_{ij}(t') \mathbf{1}_{\{i,j\} \in \mathcal{E}_{t'}} \le \tilde{y}_{ij}(t') \mathbf{1}_{\{i,j\} \in \mathcal{E}_{\delta}(t)}. $$
Then $|I_4^m|$ can be bounded as

$$
\begin{align*}
|I_4^m| &= \left| \sum_{j:j \neq m} (\hat{\pi}_j^{(m)}(t) - \pi_j^*(t)) \frac{1}{L d_\delta(t) |\mathcal{N}_{jm,\delta}(t)|} \sum_{t' \in \mathcal{N}_{jm,\delta}(t)} \sum_{l=1}^L y_{jm}^{(l)}(t') \right| \\
&\le \sum_{j:j \neq m} |\hat{\pi}_j^{(m)}(t) - \pi_j^*(t)| \left[ \frac{1}{L d_\delta(t) |\mathcal{N}_{jm,\delta}(t)|} \sum_{t' \in \mathcal{N}_{jm,\delta}(t)} \sum_{l=1}^L \tilde{y}_{jm}^{(l)}(t') \right] \mathbf{1}_{\{j,m\} \in \mathcal{E}_\delta(t)} \\
&\le \frac{N_{\max,\delta}(t)}{N_{\min,\delta}(t)} \sum_{j:j \neq m} \frac{1}{d_\delta(t)} |\hat{\pi}_j^{(m)}(t) - \pi_j^*(t)| \mathbf{1}_{\{j,m\} \in \mathcal{E}_\delta(t)}.
\end{align*}
$$

When $G_{t'} \sim \mathcal{G}(n, p(t'))$ for all $t' \in \mathcal{N}_\delta(t)$, then we know from Lemma 11 that

$$ \frac{N_{\max,\delta}(t)}{N_{\min,\delta}(t)} \le 4 \qquad (D.10) $$

holds with probability at least $1 - O(n^{-10})$. Hence, to bound $|I_4^m|$, it suffices to bound

$$ I_5^m := \sum_{j:j \neq m} \frac{1}{d_\delta(t)} |\hat{\pi}_j^{(m)}(t) - \pi_j^*(t)| \mathbf{1}_{\{j,m\} \in \mathcal{E}_\delta(t)}. $$

To this end, let us first denote $G_\delta^{(m)}(t)$ to be the union graph $G_\delta(t)$ where the item $m$ has been removed, and $\tilde{y}(t)$ to be the variables computed in (D.9). Then we have by triangle inequality that

$$ I_5^m \le \mathbb{E}[I_5^m | G_\delta^{(m)}(t), \tilde{y}] + \left| I_5^m - \mathbb{E}[I_5^m | G_\delta^{(m)}(t), \tilde{y}] \right|. \qquad (D.11) $$

The expectation can be bound using Cauchy-Schwarz inequality.

$$
\begin{align*}
\mathbb{E}[I_5^m | G_\delta^{(m)}(t), \tilde{y}] &= \sum_{j \neq m} \frac{1}{d_\delta(t)} |\hat{\pi}_j^{(m)}(t) - \pi_j^*(t)| p_\delta(t) \\
&\le \frac{\sqrt{n} p_\delta(t)}{d_\delta(t)} \|\hat{\pi}^{(m)}(t) - \pi^*(t)\|_2 \\
&\le \frac{1}{3\sqrt{n}} \|\hat{\pi}^{(m)}(t) - \pi^*(t)\|_2 \\
&\le \frac{1}{3\sqrt{n}} \|\hat{\pi}^{(m)}(t) - \hat{\pi}(t)\|_2 + \frac{1}{3\sqrt{n}} \|\hat{\pi}(t) - \pi^*(t)\|_2.
\end{align*}
$$

A bound on the first term is given by Lemma 9 and by Theorem 2 for the second term. Hence there exist constants $C_3, \tilde{C}_2 \ge 1$ so that with probability at least $1 - O(n^{-9})$, we have for all $m \in [n]$ that

$$
\begin{align*}
\mathbb{E}[I_5^m | G_\delta^{(m)}(t), \tilde{y}] & \\
&\le 64b^{\frac{5}{2}}(t) \left( \frac{4M\delta}{T\sqrt{n}} + C_3 \sqrt{\frac{\log n}{Ln^2p_\delta(t)p_{\delta,\text{sum}}(t)}} \right) \|\pi^*(t)\|_\infty + \frac{1}{3\sqrt{n}} \|\hat{\pi}(t) - \pi^*(t)\|_\infty \\
&\quad + \frac{1}{\sqrt{n}} \left( 521 \frac{M\delta nb^{\frac{7}{2}}(t)}{T} + 64\tilde{C}_2 b^{\frac{9}{2}}(t) \sqrt{\frac{1}{Ln p_\delta(t)p_{\delta,\text{sum}}(t)}} \right) \|\pi^*(t)\|_2 \\
&\le 64b^{\frac{5}{2}}(t) \left( \frac{4M\delta}{T\sqrt{n}} + C_3 \sqrt{\frac{\log n}{Ln^2p_\delta(t)p_{\delta,\text{sum}}(t)}} \right) \|\pi^*(t)\|_\infty + \frac{1}{3\sqrt{n}} \|\hat{\pi}(t) - \pi^*(t)\|_\infty
\end{align*}
$$
$$
\begin{align}
&+ \left( 512 \frac{M \delta n b^{\frac{7}{2}}(t)}{T} + 64 \tilde{C}_2 b^{\frac{9}{2}}(t) \sqrt{\frac{1}{L np_{\delta}(t) p_{\delta, \text{sum}}(t)}} \right) \|\pi^*(t)\|_{\infty} \nonumber \\
&\le \left( C'_3 \frac{M n \delta b^{\frac{7}{2}}(t)}{T} + C'_4 \frac{b^{\frac{9}{2}}(t)}{\sqrt{L np_{\delta}(t) p_{\delta, \text{sum}}(t)}} \right) \|\pi^*(t)\|_{\infty} + \frac{1}{3\sqrt{n}} \|\hat{\pi}(t) - \pi^*(t)\|_{\infty} \tag{D.12}
\end{align}
$$

for some constants $C'_3, C'_4 \ge 1$.

To bound the deviation term in (D.11), we will use Bernstein's inequality. Then, there exists a constant $c \ge 1$ such that for every $m \in [n]$,

$$
\mathbb{P} \left( \left| I_5^m - \mathbb{E}[I_5^m | G_\delta^{(m)}(t), \tilde{y}] \right| \ge c \frac{\sqrt{np_\delta(t) \log n} + \log n}{d_\delta(t)} \|\hat{\pi}^{(m)}(t) - \pi^*(t)\|_{\infty} \right) \le 2n^{-10}.
$$

Hence using a union bound over $[n]$, we see that there exists a suitably large constant $c' \ge 1$ so
that with probability at least $1 - O(n^{-9})$,

$$
\left| I_5^m - \mathbb{E}[I_5^m | G_\delta^{(m)}(t), \tilde{y}] \right| \le c' \sqrt{\frac{\log n}{np_\delta(t)}} \|\hat{\pi}^{(m)}(t) - \pi^*(t)\|_{\infty} \quad \forall m \in [n]. \quad (D.13)
$$

It remains to bound $\|\hat{\pi}^{(m)}(t) - \pi^*(t)\|_{\infty}$. Using Lemma 9 along with the triangular inequality

$$
\| \hat{\pi}^{(m)}(t) - \pi^{*}(t) \|_{\infty} \le \| \hat{\pi}^{(m)}(t) - \pi(t) \|_{2} + \| \hat{\pi}(t) - \pi^{*}(t) \|_{\infty},
$$

we have that with probability at least $1 - O(n^{-9})$,

$$
\|\hat{\pi}^{(m)}(t) - \pi^{*}(t)\|_{\infty} \le 192b^{\frac{5}{2}}(t) \left( \frac{4M\delta}{T} + C_3 \sqrt{\frac{\log n}{Lnp_{\delta}(t)p_{\delta,\text{sum}}(t)}} \right) \|\pi^{*}(t)\|_{\infty} + 2\|\hat{\pi}(t) - \pi^{*}(t)\|_{\infty}.
$$

Upon plugging this in (D.13), the latter simplifies to

$$
\begin{align}
& \left| I_5^m - \mathbb{E}[I_5^m | G_\delta^{(m)}(t), \tilde{y}] \right| \nonumber \\
& \le C_3'' \sqrt{\frac{\log n}{np_\delta(t)}} \|\hat{\pi}(t) - \pi^*(t)\|_{\infty} + C_4'' b^{\frac{5}{2}}(t) \sqrt{\frac{\log n}{np_\delta(t)}} \left( \frac{4M\delta}{T} + \sqrt{\frac{\log n}{Lnp_\delta(t)p_{\delta,\text{sum}}(t)}} \right) \|\pi^*(t)\|_{\infty} \tag{D.14}
\end{align}
$$

for some constants $C_3'', C_4'' \ge 1$.

Finally, combining (D.10), (D.11), (D.12) and (D.14), we conclude that with probability at least
$1 - O(n^{-9})$, it holds for all $m \in [n]$ that

$$
|I_4^m| \le \left( C_7 \frac{M n \delta b^{\frac{7}{2}}(t)}{T} + C_8 \frac{b^{\frac{5}{2}}(t) \max\left\{b^2(t), \frac{\log n}{\sqrt{np_\delta(t)p_{\delta,\text{sum}}(t)}}\right\}}{\sqrt{L np_\delta(t) p_{\delta,\text{sum}}(t)}} \right) \|\pi^*(t)\|_{\infty} + C_9 \sqrt{\frac{\log n}{np_\delta(t)}} \|\hat{\pi}(t) - \pi^*(t)\|_{\infty}
$$

for some constants $C_7, C_8, C_9 \ge 1$.