Spectral Method and Regularized MLE Are Both Optimal for Top-$K$ Ranking

Introduction

Imagine we have a large collection of $n$ items, and we are given partially revealed comparisons between pairs of items. These paired comparisons are collected in a non-adaptive fashion, and could be highly noisy and incomplete. The aim is to aggregate these partial preferences so as to identify the $K$ items that receive the highest ranks. This problem, which is called top-$K$ rank aggregation, finds applications in numerous contexts, including web search (Dwork et al., 2001), recommendation systems (Baltrunas et al., 2010), sports competition (Masse, 1997), to name just a few. The challenge is both statistical and computational: how can one achieve reliable top-$K$ ranking from a minimal number of pairwise comparisons, while retaining computational efficiency?

To address the aforementioned challenge, many prior approaches have been put forward based on certain statistical models. Arguably one of the most widely used parametric models is the Bradley-Terry-Luce (BTL) model (Bradley and Terry, 1952; Luce, 1959), which assigns a latent preference score $\{w_{i}^{*}\}_{1\leq i\leq n}$ to each of the $n$ items. The BTL model posits that: the chance of each item winning a paired comparison is determined by the relative scores of the two items involved, or more precisely,

in each comparison of item $i$ against item $j$. The items are repeatedly compared in pairs according to this parametric model. The task then boils down to identifying the $K$ items with the highest preference scores, given these pairwise comparisons.

Among the ranking algorithms tailored to the BTL model, the following two procedures have received particular attention, both of which rank the items based on appropriate estimates of the latent preference scores.

The spectral method. By connecting the winning probability in (1) with the transition probability of a reversible Markov chain, the spectral method attempts recovery of $\left\{w_{i}^{*}\right\}$ via the leading left eigenvector of a sample transition matrix. This procedure, also known as Rank Centrality (Negahban et al., 2017a), bears similarity to the PageRank algorithm.

The maximum likelihood estimator (MLE). This approach proceeds by finding the score assignment that maximizes the likelihood function (Ford, 1957). When parameterized appropriately, solving the MLE becomes a convex program, and hence is computationally feasible. There are also important variants of the MLE that enforce additional regularization.

As we will elaborate later, the spectral method part of the preceding questions was recently explored by (Jang et al., 2016), for a regime where a relatively large fraction of item pairs have been compared. However, it remains unclear how well the spectral method can perform in a much broader — and often much more challenging — regime, where the fraction of item pairs being compared may be vanishingly small. Additionally, the ranking accuracy of the MLE (and its variants) remains unknown.

2 Main contributions

The central focal point of the current paper is to assess the accuracy of both the spectral method and the regularized MLE in top-$K$ identification. Assuming that the pairs of items being compared are randomly selected and that the preference scores fall within a fixed dynamic range, our paper delivers a somewhat surprising message:

Both the spectral method and the regularized MLE achieve perfect identification of top-$K$ ranked items under optimal sample complexity (up to some constant factor)!

3 Notation

Additionally, the notation $f(n)=O\left(g(n)\right)$ or $f(n)\lesssim g(n)$ means there is a constant $c>0$ such that $\left|f(n)\right|\leq c|g(n)|$, $f(n)=\Omega\left(g(n)\right)$ or $f(n)\gtrsim g(n)$ means there is a constant $c>0$ such that $|f(n)|\geq c\left|g(n)\right|$, $f(n)=\Theta\left(g(n)\right)$ or $f(n)\asymp g(n)$ means that there exist constants $c_{1},c_{2}>0$ such that $c_{1}|g(n)|\leq|f(n)|\leq c_{2}|g(n)|$, and $f(n)=o(g(n))$ means $\lim_{n\rightarrow\infty}\frac{f(n)}{g(n)}=0$.

Statistical models and main results

We begin with a formal introduction of the Bradley-Terry-Luce parametric model for binary comparisons.

Preference scores. As introduced earlier, we assume the existence of a positive latent score vector

that comprises the underlying preference scores $\left\{w_{i}^{*}>0\right\}_{1\leq i\leq n}$ assigned to each of the $n$ items. Alternatively, it is sometimes more convenient to reparameterize the score vector by

These scores are assumed to fall within a dynamic range given by

for all $1\leq i\leq n$ and for some $w_{\min}>0,$ $w_{\max}>0$, $\theta_{\min}=\log w_{\min}$, and $\theta_{\max}=\log w_{\max}$. We also introduce the condition number as

Notably, the current paper primarily focuses on the case with a fixed dynamic range (i.e. $\kappa$ is a fixed constant independent of $n$), although we will also discuss extensions to the large dynamic range regime in Section 3. Without loss of generality, it is assumed that

meaning that items $1$ through $K$ are the desired top-$K$ ranked items.

Comparison graph. Let $\mathcal{G}=(\mathcal{V},\mathcal{E})$ stand for a comparison graph, where the vertex set $\mathcal{V}=\left\{1,2,\ldots,n\right\}$ represents the $n$ items of interest. The items $i$ and $j$ are compared if and only if $(i,j)$ falls within the edge set $\mathcal{E}$. Unless otherwise noted, we assume that $\mathcal{G}$ is drawn from the Erdős–Rényi random graph $\mathcal{G}_{n,p}$, such that an edge between any pair of vertices is present independently with some probability $p$. In words, $p$ captures the fraction of item pairs being compared.

By convention, we set $y_{i,j}^{(l)}=1-y_{j,i}^{(l)}$ for all $(i,j)\in\mathcal{E}$ throughout the paper. This is also known as the logistic pairwise comparison model, due to its strong resemblance to logistic regression. It is self-evident that the sufficient statistics under this model are given by

To simplify the notation, we shall also take

Goal. The goal is to identify the set of top-$K$ ranked items — that is, the set of $K$ items that enjoy the largest preference scores — from the pairwise comparison data $\bm{y}$.

2 Algorithms

The spectral ranking algorithm, or Rank Centrality (Negahban et al., 2017a), is motivated by the connection between the pairwise comparisons and a random walk over a directed graph. The algorithm starts by converting the pairwise comparison data $\bm{y}$ into a transition matrix $\bm{P}=[P_{i,j}]_{1\leq i,j\leq n}$ in such a way that

To develop some intuition regarding why this spectral algorithm gives a reasonable estimate of $\bm{w}^{*}$, it is perhaps more convenient to look at the population transition matrix $\bm{P}^{*}=[P_{i,j}^{*}]_{1\leq i,j\leq n}$:

which coincides with $\bm{P}$ by taking $L\rightarrow\infty$. It can be seen that the normalized score vector

is the stationary distribution of the Markov chain induced by the transition matrix $\bm{P}^{*}$, since $\bm{P}^{*}$ and $\bm{\pi}^{*}$ are in detailed balance, namely,

As a result, one expects the stationary distribution of the sample version $\bm{P}$ to form a good estimate of $\bm{w}^{*}$, provided the sample size is sufficiently large.

2.2 The regularized MLE

Under the BTL model, the negative log-likelihood function conditioned on $\mathcal{G}$ is given by (up to some global scaling)

The regularized MLE then amounts to solving the following convex program

for a regularization parameter $\lambda>0$. As will be discussed later, we shall adopt the choice $\lambda\asymp\sqrt{\frac{np\log n}{L}}$ throughout this paper. For the sake of brevity, we let $\bm{\theta}$ represent the resulting penalized maximum likelihood estimate whenever it is clear from the context. Similar to the spectral method, one reports the $K$ items associated with the $K$ largest entries of $\bm{\theta}$.

3 Main results

The most challenging part of top-$K$ ranking is to distinguish the $K$-th and the $(K+1)$-th items. In fact, the score difference of these two items captures the distance between the item sets $\left\{1,\ldots,K\right\}$ and $\left\{K+1,\ldots,n\right\}$. Unless their latent scores are sufficiently separated, the finite-sample nature of the model would make it infeasible to distinguish these two critical items. With this consideration in mind, we define the following separation measure

This metric turns out to play a crucial role in determining the minimal sample complexity for perfect top-$K$ identification.

The main finding of this paper concerns the optimality of both the spectral method and the regularized MLE in the presence of a fixed dynamic range (i.e. $\kappa=O(1)$). Recall that under the BTL model, the total number $N$ of samples we collect concentrates sharply around its mean, namely,

occurs with high probability. Our main result is stated in terms of the sample complexity required for exact top-$K$ identification.

Consider the pairwise comparison model specified in Section 2.1 with $\kappa=O(1)$. Suppose that $p>\frac{c_{0}\log n}{n}$ and that

for some sufficiently large positive constants $c_{0}$ and $c_{1}$. Further assume $L\leq c_{2}\cdot n^{c_{3}}$ for any absolute constants $c_{2},c_{3}>0$. With probability exceeding $1-O(n^{-5})$, the set of top-$K$ ranked items can be recovered exactly by the spectral method given in Algorithm 1, and by the regularized MLE given in (13). Here, we take $d=c_{d}np$ in the spectral method and $\lambda=c_{\lambda}\sqrt{\frac{np\log n}{L}}$ in the regularized MLE, where $c_{d}\geq 2$ and $c_{\lambda}>0$ are some absolute constants.

We emphasize that $p\geq\frac{c_{0}\log n}{n}$ for $c_{0}\geq 1$ is a fundamental requirement for the ranking task. In fact, if $p<(1-\epsilon)\frac{\log n}{n}$ for any constant $\epsilon>0$, then the comparison graph $\mathcal{G}\sim\mathcal{G}_{n,p}$ is disconnected with high probability. This means that there exists at least one isolated item (which has not been compared with any other item) and cannot be ranked.

In fact, the assumption that $L\leq c_{2}\cdot n^{c_{3}}$ for any absolute constants $c_{2},c_{3}>0$ is not needed for the spectral method.

Here, we assume the same number of comparisons $L$ to simplify the presentation as well as the proof. The result still holds true if we have distinct $L_{i,j}$’s for each $i\neq j$, as long as $n^{2}p\min_{i\neq j}L_{i,j}\gtrsim\frac{n\log n}{\Delta_{K}^{2}}$.

Theorem 1 asserts that both the spectral method and the regularized MLE achieve a sample complexity on the order of $\frac{n\log n}{\Delta_{{\it K}}^{2}}$. Encouragingly, this sample complexity coincides with the minimax limit identified in (Chen and Suh, 2015, Theorem 2) in the fixed dynamic range, i.e. $\kappa=O(1)$.

Fix $\epsilon\in(0,\frac{1}{2})$, and suppose that

where $c_{2}=w_{\min}^{4}/(4w_{\max}^{4})$. Then for any ranking procedure $\psi$, one can find a score vector $\bm{w}^{*}$ with separation $\Delta_{K}$ such that $\psi$ fails to retrieve the top-$K$ items with probability at least $\epsilon$.

We are now positioned to compare our results with Jang et al. (2016), which also investigates the accuracy of the spectral method for top-$K$ ranking. Specifically, Theorem 3 in Jang et al. (2016) establishes the optimality of the spectral method for the relatively dense regime where

In this regime, however, the total sample size necessarily exceeds

which rules out the possibility of achieving minimal sample complexity if $\Delta_{K}$ is sufficiently large. For instance, consider the case where $\Delta_{K}\asymp 1$, then the optimal sample size — as revealed by Theorem 1 or (Chen and Suh, 2015, Theorem 1) — is on the order of

which is a factor of $\sqrt{\frac{n}{\log n}}$ lower than the bound in (18). By contrast, our results hold all the way down to the sparsest possible regime where $p\asymp\frac{\log n}{n}$, confirming the optimality of the spectral method even for the most challenging scenario. Furthermore, we establish that the regularized MLE shares the same optimality guarantee as the spectral method, which was previously out of reach.

4 Optimal control of entrywise estimation errors

Consider the pairwise comparison model in Section 2.1 with $\kappa=O(1)$. Suppose $p>\frac{c_{0}\log n}{n}$ for some sufficiently large constant $c_{0}>0$. Choose $d=c_{d}np$ for some constant $c_{d}\geq 2$ in Algorithm 1. Then the spectral estimate $\bm{\pi}$ satisfies

with probability $1-O(n^{-5})$, where $\bm{\pi}^{*}$ is the normalized score vector (cf. (10)).

Consider the pairwise comparison model specified in Section 2.1 with $\kappa=O\left(1\right)$. Suppose that $p\geq\frac{c_{0}\log n}{n}$ for some sufficiently large constant $c_{0}>0$ and that $L\leq c_{2}\cdot n^{c_{3}}$ for any absolute constants $c_{2},c_{3}>0$. Set the regularization parameter to be $\lambda=c_{\lambda}\sqrt{\frac{np\log n}{L}}$ for some absolute constant $c_{\lambda}>0$. Then the regularized MLE $\bm{\theta}$ satisfies

with probability exceeding $1-O(n^{-5})$, where $\overline{\theta}^{*}:=\frac{1}{n}\bm{1}^{\top}\bm{\theta}^{*}$ and $e^{\bm{\theta}}:=[e^{\theta_{1}},\cdots,e^{\theta_{n}}]^{\top}$.

with high probability. Similar theoretical guarantees have been derived for another variant of the MLE (the constrained version) under a uniform sampling model as well (Negahban et al., 2017a). In comparison, our results indicate that the estimation errors for both algorithms are almost evenly spread out across all coordinates rather than being localized or clustered. Notably, the pointwise errors revealed by Theorems 3-4 immediately lead to exact top-$K$ identification as claimed by Theorem 1.

In what follows, we prove the theorem for the spectral method part. The regularized MLE part follows from an almost identical argument and hence is omitted.

Since the spectral algorithm ranks the items in accordance with the score estimate $\bm{\pi}$, it suffices to demonstrate that

To this end, we first apply the triangle inequality to get

In addition, it follows from Theorem 3 as well as our sample complexity assumption that

These conditions taken collectively imply that $\frac{\|\bm{\pi}-\bm{\pi}^{*}\|_{\infty}}{\|\bm{\pi}^{*}\|_{\infty}}<\frac{1}{2}\Delta_{K}$ as long as $\frac{npL\Delta_{K}^{2}}{\log n}$ exceeds some sufficiently large constant. Substitution into (21) reveals that $\pi_{i}-\pi_{j}>0$, as claimed. ∎

5 Heuristic arguments

We pause to develop some heuristic explanation as to why the estimation errors are expected to be spread out across all entries. For simplicity, we focus on the case where $p=1$ and $L$ is sufficiently large, so that $\bm{y}$ and $\bm{P}$ sharply concentrate around $\bm{y}^{*}$ and $\bm{P}^{*}$, respectively.

We begin with the spectral algorithm. Since $\bm{\pi}$ and $\bm{\pi}^{*}$ are respectively the invariant distributions of the Markov chains induced by $\bm{P}$ and $\bm{P}^{*}$, we can decompose

When $p=1$ and $\frac{w_{\max}}{w_{\min}}\asymp 1$, the entries of $\bm{\pi}^{*}$ (resp. the off-diagonal entries of $\bm{P}^{*}$ and $\bm{P}-\bm{P}^{*}$) are all of the same order and, as a result, the energy of the uncertainty term $\bm{\xi}$ is spread out (using standard concentration inequalities). In fact, we will demonstrate in Section 5.2 that

which coincides with the optimal rate. Further, if we look at each entry of (22), then for all $1\leq m\leq n$,

By construction of the transition matrix, one can easily verify that $P_{m,m}$ is bounded away from $1$ and $P_{j,m}\asymp\frac{1}{n}$ for all $j\neq m$. As a consequence, the identity $\bm{\pi}^{\top}\bm{P}=\bm{\pi}^{\top}$ allows one to treat each $\pi_{m}-\pi_{m}^{*}$ as a mixture of three effects: (i) the first term of (30) behaves as an entrywise contraction of the error; (ii) the second term of (30) is a (nearly uniformly weighted) average of the errors over all coordinates, which can essentially be treated as a smoothing operator applied to the error components; and (iii) the uncertainty term $\xi_{m}$. Rearranging terms in (30), we are left with

There are two possibilities compatible with this bound (32): (1) $\left\|\bm{\pi}-\bm{\pi}^{*}\right\|_{\infty}\lesssim\frac{1}{n}\sum_{i=1}^{n}\left|\pi_{i}-\pi_{i}^{*}\right|$, and (2) $\left\|\bm{\pi}-\bm{\pi}^{*}\right\|_{\infty}\lesssim\|\bm{\xi}\|_{\infty}\lesssim\frac{\|\bm{\pi}-\bm{\pi}^{*}\|_{2}}{\|\bm{\pi}^{*}\|_{2}}\left\|\bm{\pi}^{*}\right\|_{\infty}$ by (23). In either case, the errors are fairly delocalized, revealing that

We now move on to the regularized MLE, following a very similar argument. By the optimality condition that $\nabla\mathcal{L}_{\lambda}\left(\bm{\theta}\right)=\bm{0}$, one can derive (for some $\eta$ to be specified later)

Write $\nabla^{2}\mathcal{L}_{\lambda}\left(\bm{\theta}^{*}\right)=\bm{D}-\bm{A}$, where $\bm{D}$ and $\bm{A}$ denote respectively the diagonal and off-diagonal parts of $\nabla^{2}\mathcal{L}_{\lambda}\left(\bm{\theta}^{*}\right)$. Under our assumptions, one can check that $D_{m,m}\asymp n$ for all $1\leq m\leq n$ and $A_{j,m}\asymp 1$ for any $j\neq m$. With these notations in place, one can write the entrywise error as follows

By choosing $\eta=c_{2}/n$ for some sufficiently small constant $c_{2}>0$ , we get $1-\eta D_{m,m}<1$ and $\eta A_{j,m}\asymp 1/n$. Therefore, the right-hand side of the above relation also comprises a contraction term as well as an error smoothing term, similar to (30). Carrying out the same argument as for the spectral method, we see that the estimation errors of the regularized MLE are expected to be spread out.

6 Numerical experiments

It is worth noting that extensive numerical experiments on both synthetic and real data have already been conducted in Negahban et al. (2017a) to confirm the practicability of both the spectral method and the regularized MLE. See also Chen and Suh (2015) for the experiments on the Spectral-MLE algorithm. This section provides some additional simulations to complement their experimental results as well as our theory. Throughout the experiments, we set the number of items $n$ to be $200$, while the number of repeated comparisons $L$ and the edge probability $p$ can vary with the experiments. Regarding the tuning parameters, we choose $d=2d_{\max}$ in the spectral method where $d_{\max}$ is the maximum degree of the graph and $\lambda=2\sqrt{\frac{np\log n}{L}}$ in the regularized MLE, which are consistent with the configurations considered in the main theorems. Additionally, we also display the experimental results for the unregularized MLE, i.e. $\lambda=0$. All of the results are averaged over 100 Monte Carlo simulations.

Further, we examine the top-$K$ ranking accuracy of all three methods. Here, we fix $p=0.25$ and $L=20$, set ${\it K}=10$, and let $w^{*}_{i}=1$ for all $1\leq i\leq K$ and $w^{*}_{j}=1-\Delta$ for all $K+1\leq j\leq n$. By construction, the score separation satisfies $\Delta_{K}=\Delta$. Figure 3 illustrates the accuracy in identifying the top-$K$ ranked items. The performance of them improves when the score separation becomes larger, which matches our theory in Theorem 1.

7 Other related works

The problem of ranking based on partial preferences has received much attention during the past decade. Two types of observation models have been considered: $(1)$ the cardinal-based model, where users provide explicit numerical ratings of the items, $(2)$ the ordinal-based model, where users are asked to make comparative measurements. See Ammar and Shah (2011) for detailed comparisons between them.

In terms of the ordinal-based model — and in particular, ranking from pairwise comparisons — both parametric and nonparametric models have been extensively studied. For example, Hunter (2004) examined variants of the parametric BTL model, and established the convergence properties of the minorization-maximization algorithm for computing the MLE. Moreover, the BTL model falls under the category of low-rank parametric models, since the preference matrix is generated by passing a rank-2 matrix through the logistic link function (Rajkumar and Agarwal, 2016). Additionally, the work Jiang et al. (2011) proposed a least-squares type method to estimate the full ranking, which generalizes the simple Borda count algorithm (Ammar and Shah, 2011). For many of these algorithms, the sample complexities needed for perfect total ranking were determined by Rajkumar and Agarwal (2014), although the top-$K$ ranking accuracy was not considered there.

Going beyond the parametric models, a recent line of works Shah et al. (2017); Shah and Wainwright (2015); Chen et al. (2017); Pananjady et al. (2017) considered the nonparametric stochastically transitive model, where the only model assumption is that the comparison probability matrix follows certain transitivity rules. This type of models subsumes the BTL model as a special case. For instance, Shah and Wainwright (2015) suggested a simple counting-based algorithm which can reliably recover the top-$K$ ranked items for various models. However, the sampling paradigm considered therein is quite different from ours in the sparse regime; for instance, their model does not come close to the setting where $p$ is small but $L$ is large, which is the most challenging regime of the model adopted in our paper and Negahban et al. (2017a); Chen and Suh (2015).

All of the aforementioned papers concentrate on the case where there is a single ground-truth ordering. It would also be interesting to investigate the scenarios where different users might have different preference scores. To this end, Negahban et al. (2017b); Lu and Negahban (2014) imposed the low-rank structure on the underlying preference matrix and adopted the nuclear-norm relaxation approach to recover the users’ preferences. Additionally, several papers explored the ranking problem for the more general Plackett-Luce model (Hajek et al., 2014; Soufiani et al., 2013), in the presence of adaptive sampling (Jamieson and Nowak, 2011; Busa-Fekete et al., 2013; Heckel et al., 2016; Agarwal et al., 2017), for the crowdsourcing scenario (Chen et al., 2013), and in the adversarial setting (Suh et al., 2017). These are beyond the scope of the present paper.

Finally, the family of spectral methods has been successfully applied in numerous applications, e.g. matrix completion (Keshavan et al., 2010), phase retrieval (Chen and Candès, 2017), graph clustering (Rohe et al., 2011; Abbe et al., 2017), joint alignment (Chen and Candes, 2016). All of them are designed based on the eigenvectors of some symmetric matrix, or the singular vectors if the matrix of interest is asymmetric. Our paper contributes to this growing literature by establishing a sharp eigenvector perturbation analysis framework for an important class of asymmetric matrices — the probability transition matrices.

Extension: general dynamic range

All of the preceding results concern the regime with a fixed dynamic range (i.e. $\kappa=O(1)$). This section moves on to discussing the case with large $\kappa$.

To start with, by going through the same proof technique, we can readily obtain — in the general $\kappa$ setting — the following performance guarantees for both the spectral estimate $\bm{\pi}$ and the regularized MLE $\bm{\theta}$.

Consider the pairwise comparison model in Section 2.1. Suppose that $p>\frac{c_{0}\kappa^{5}\log n}{n}$ for some sufficiently large constant $c_{0}>0$, and choose $d=c_{d}np$ for some constant $c_{d}\geq 2$ in Algorithm 1. Then with probability exceeding 1-O\big{(}n^{-5}\big{)},

the spectral estimate $\bm{\pi}$ satisfies

where $\bm{\pi}^{*}$ is the normalized score vector as defined in (10).

the set of top-$K$ ranked items can be recovered exactly by the spectral method given in Algorithm 1, as long as

for some sufficiently large constant $c_{1}>0$.

Consider the pairwise comparison model in Section 2.1. Suppose that $p\geq\frac{c_{0}\kappa^{4}\log n}{n}$ for some sufficiently large constant $c_{0}>0$ and that $L\leq c_{2}\cdot n^{c_{3}}$ for any absolute constants $c_{2},c_{3}>0$. Set the regularization parameter to be $\lambda=c_{\lambda}\frac{1}{\log\kappa}\sqrt{\frac{np\log n}{L}}$ for some absolute constant $c_{\lambda}>0$. Then with probability exceeding 1-O\big{(}n^{-5}\big{)},

the regularized MLE $\bm{\theta}$ satisfies

where $\overline{\theta}^{*}:=\frac{1}{n}\bm{1}^{\top}\bm{\theta}^{*}$ and $e^{\bm{\theta}}:=[e^{\theta_{1}},\cdots,e^{\theta_{n}}]^{\top}$.

the set of top-$K$ ranked items can be recovered exactly by the regularized MLE given in (13), as long as

for some sufficiently large constant $c_{1}>0$.

Notably, the achievability bounds for top-$K$ ranking in Theorems 5–6 do not match the lower bound asserted in Theorem 2 in terms of $\kappa$. This is partly because the separation measure $\Delta_{K}$ fails to capture the information bottleneck for the general $\kappa$ setting. In light of this, we introduce the following new measure that seems to be a more suitable metric to reflect the hardness of the top-$K$ ranking problem:

which will be termed the generalized separation measure. Informally, $(\Delta_{K}^{*})^{2}$ is a reasonably tight upper bound on certain normalized KL divergence metric. With this metric in place, we derive another lower bound as follows.

Fix $\epsilon\in(0,\frac{1}{2})$, and let $\mathcal{G}\sim\mathcal{G}_{n,p}$. Consider any preference score vector $\bm{w}^{*}$, and let $\Delta_{K}^{*}$ denote its generalized separation. If

The preceding sample complexity lower bound scales inversely proportionally to $\left(\Delta_{K}^{*}\right)^{2}$. To see why this generalized measure may be more suitable compared to the original separation metric, we single out three examples in Appendix B. Unfortunately, our current analyses do not yield a matching upper bound with respect to $\Delta_{K}^{*}$ unless $\kappa$ is a constant. For instance, the analysis of the spectral method relies on the eigenvector perturbation bound (Theorem 8), where the spectral gap and matrix perturbation play a crucial rule. However, the current results for controlling these quantities have explicit dependency on $\kappa$ Negahban et al. (2017a). It is not clear whether we could incorporate the new measure to eliminate such dependency on $\kappa$. This calls for more refined analysis techniques, which we leave for future investigation.

Moreover, it is not obvious whether the spectral method alone or the regularized MLE alone can achieve the minimal sample complexity in the general $\kappa$ regime. It is possible that one needs to first screen out those items with extremely high or low scores using methods like Borda count (Ammar and Shah, 2012), as advocated by (Negahban et al., 2017a; Chen and Suh, 2015; Jang et al., 2016). All in all, finding tight upper bounds for general $\kappa$ remains an open question.

Discussion

This paper justifies the optimality of both the spectral method and the regularized MLE for top-$K$ rank aggregation for the fixed dynamic range case. Our theoretical studies are by no means exhaustive, and there are numerous directions that would be of interest for future investigations. We point out a few possibilities as follows.

General condition number $\kappa$. As mentioned before, our current theory is optimal in the presence of a fixed dynamic range with $\kappa=O\left(1\right)$. We have also made a first attempt in considering the large $\kappa$ regime. It is desirable to characterize the statistical and computational limits for more general $\kappa$.

Goodness-of-fit. Throughout this paper, we have assumed the BTL model captures the randomness underlying the data we collect. A practical question is whether the real data actually follows the BTL model. It would be interesting to investigate how to test the goodness-of-fit of this model.

Unregularized MLE. We have studied the optimality of the regularized MLE with the regularization parameter $\lambda\asymp\sqrt{\frac{np\log n}{L}}$. Our analysis relies on the regularization term to obtain convergence of the gradient descent algorithm (see Lemma 11). It is natural to ask whether such a regularization term is necessary or not. This question remains open.

More general comparison graphs. So far we have focused on a tractable but somewhat restrictive comparison graph, namely, the Erdős–Rényi random graph. It would certainly be important to understand the performance of both methods under a broader family of comparison graphs, and to see which algorithms would enable optimal sample complexities under general sampling patterns.

Analysis for the spectral method

This section is devoted to proving Theorem 5 and hence Theorem 3, which characterizes the pointwise error of the spectral estimate.

Here, we gather some preliminary facts about reversible Markov chains as well as the Erdős–Rényi random graph.

The first important result concerns the eigenvector perturbation for probability transition matrices, which can be treated as the analogue of the celebrated Davis-Kahan $\sin\Theta$ theorem (Davis and Kahan, 1970). Due to its potential importance for other problems, we promote it to a theorem as follows.

Suppose that $\bm{P}$, $\hat{\bm{P}}$, and $\bm{P}^{*}$ are probability transition matrices with stationary distributions $\bm{\pi}$, $\hat{\bm{\pi}}$, $\bm{\pi}^{*}$, respectively. Also, assume that $\bm{P}^{*}$ represents a reversible Markov chain. When \big{\|}\bm{P}-\hat{\bm{P}}\big{\|}_{\bm{\pi}^{*}}<1-\max\left\{\lambda_{2}(\bm{P}^{*}),-\lambda_{n}\left(\bm{P}^{*}\right)\right\}, it holds that

Consider the pairwise comparison model specified in Section 2.1 with $\kappa=O(1)$. Suppose $p\geq c_{0}\frac{\log n}{n}$ for some sufficiently large constant $c_{0}>0$ and $d\geq c_{d}np$ for $c_{d}\geq 2$ in Algorithm 1. With probability exceeding $1-O(n^{-5})$, one has

The next result is concerned with the concentration of the vertex degrees in an Erdős–Rényi random graph.

Suppose that $\mathcal{G}\sim\mathcal{G}_{n,p}$. Let $d_{i}$ be the degree of node $i$, $d_{\min}=\min_{1\leq i\leq n}d_{i}$ and $d_{\max}=\max_{1\leq i\leq n}d_{i}$. If $p\geq\frac{c_{0}\log n}{n}$ for some sufficiently large constant $c_{0}>0$, then the following event

The proof follows from the standard Chernoff bound and is hence omitted. ∎

Since $d$ is chosen to be $c_{d}np$ for some constant $c_{d}\geq 2$, we have, by Lemma 1, that the maximum vertex degree obeys $d_{\max}<d$ with high probability.

2 Proof outline of Theorem 5

In this subsection, we outline the proof of Theorem 5.

Recall that $\bm{\pi}=\left[\pi_{1},\cdots,\pi_{n}\right]^{\top}$ and $\bm{\pi}^{*}=\left[\pi_{1}^{*},\cdots,\pi_{n}^{*}\right]^{\top}$ are the stationary distributions associated with $\bm{P}$ and $\bm{P}^{*}$, respectively. This gives

For each $1\leq m\leq n$, one can decompose

where $\bm{P}_{\cdot m}$ (resp. $\bm{P}^{*}_{\cdot m}$) denotes the $m$-th column of $\bm{P}$ (resp. $\bm{P}^{*}$). Then it boils down to controlling $I_{1}^{m}$, $I_{2}^{m}$ and $\sum_{j:j\neq m}(\pi_{j}-\pi_{j}^{*})P_{j,m}$.

Since $\bm{\pi}^{*}$ is deterministic while $\bm{P}$ is random, we can easily control $I_{1}^{m}$ using Hoeffding’s inequality. The bound is the following.

With probability exceeding $1-O(n^{-5})$, one has

Next, we show the term $I_{2}^{m}$ behaves as a contraction of $|\pi_{m}-\pi_{m}^{*}|$.

With probability exceeding $1-O(n^{-5})$, there exists some constant $c>0$ such that for all $1\leq m\leq n$,

The statistical dependency between $\bm{\pi}$ and $\bm{P}$ introduces difficulty in obtaining a sharp estimate of the third term $\sum_{j:j\neq m}(\pi_{j}-\pi_{j}^{*})P_{j,m}$. Nevertheless, the leave-one-out technique helps us decouple the dependency and obtain effective control of this term. The key component of the analysis is the introduction of a new probability transition matrix $\bm{P}^{(m)}$, which is a leave-one-out version of the original matrix $\bm{P}$. More precisely, $\bm{P}^{(m)}$ replaces all of the transition probabilities involving the $m$-th item with their expected values (unconditional on $\mathcal{G}$); that is, for any $i\neq j$,

with $y_{i,j}^{*}:=\frac{w_{j}^{*}}{w_{i}^{*}+w_{j}^{*}}$. For any $1\leq i\leq n$, set

in order to ensure that $\bm{P}^{(m)}$ is a probability transition matrix. In addition, we let $\bm{\pi}^{\left(m\right)}$ be the stationary distribution of the Markov chain induced by $\bm{P}^{(m)}$. As will be demonstrated later, the main advantages of introducing $\bm{\pi}^{(m)}$ are two-fold: (1) the original spectral estimate $\bm{\pi}$ is very well approximated by $\bm{\pi}^{(m)}$, and (2) $\bm{\pi}^{(m)}$ is statistically independent of the connectivity of the $m$-th node and the comparisons with regards to the $m$-th item. Now we further decompose $\sum_{j:j\neq m}(\pi_{j}-\pi_{j}^{*})P_{j,m}$:

For $I_{3}^{m}$, we apply the Cauchy-Schwarz inequality to obtain that with probability at least $1-O(n^{-10}),$

Suppose that $np\gamma^{2}>c\kappa\log n$ for some sufficiently large constant $c>0$. With probability at least $1-O\left(n^{-5}\right)$,

where $\kappa=w_{\max}/w_{\min}$ and $\gamma=1-\max\left\{\lambda_{2}(\bm{P}^{*}),-\lambda_{n}\left(\bm{P}^{*}\right)\right\}-\|\bm{P}-\bm{P}^{*}\|_{\bm{\pi}^{*}}$.

Using (Negahban et al., 2017a, Lemmas 3 and 4) and our Lemma 1, we can bound $\gamma$ from below:

Under the model specified in Section 2.1, if $p\geq c_{0}\frac{\log n}{n}\max\{1,\frac{\kappa^{5}}{L}\}$ for some sufficiently large constant $c_{0}>0$, then with probability at least $1-O(n^{-5})$,

In order to control $I_{4}^{m}$, we exploit the statistical independence between $\bm{\pi}^{(m)}$ and $\bm{P}_{\cdot m}$. Specifically, we demonstrate that:

Suppose that $p>\frac{c_{0}\log n}{n}$ for some sufficiently large constant $c_{0}>0$. With probability at least $1-O\left(n^{-10}\right)$,

The above bound depends on both $\|\bm{\pi}^{(m)}-\bm{\pi}\|_{2}$ and $\|\bm{\pi}^{(m)}-\bm{\pi}^{*}\|_{\infty}$. We can invoke Lemma 4 and the inequality $\|\bm{\pi}^{(m)}-\bm{\pi}^{*}\|_{\infty}\leq\|\bm{\pi}^{(m)}-\bm{\pi}\|_{2}+\|\bm{\pi}-\bm{\pi}^{*}\|_{\infty}$ to reach

Finally we put the preceding bounds together. When $\frac{np}{\kappa^{5}\log n}$ is large enough, with high probability, for some absolute constants $c_{1},c_{2},c_{3}>0$ one has

simultaneously for all $1\leq m\leq n$. By taking the maximum over $m$ on the left-hand side and combining terms, we get

Hence, as long as $\frac{np}{\kappa^{5}\log n}$ is sufficiently large, one has

which further leads to $\alpha_{1}\gtrsim 1/\kappa$, $\alpha_{2}\lesssim 1$, and

This finishes the proof of Theorem 5 and Theorem 3.

Analysis for the regularized MLE

This combined with the fact that $\theta_{\max}-\theta_{\min}=\log\kappa$ reveals that

In addition, we assume that $L=O\left(n^{5}\right)$ in this section. It is straightforward to extend the proof to cover $L\leq c_{2}\cdot n^{c_{3}}$ for any constants $c_{2},c_{3}>0$.

Before proceeding to the proof, we gather some basic facts. To begin with, the gradient and the Hessian of $\mathcal{L}\left(\cdot;\bm{y}\right)$ in (12) can be computed as

Let $\lambda$ be as specified in Theorem 6. The following event

occurs with probability exceeding $1-O(n^{-10})$.

The following lemmas characterize the smoothness and the strong convexity of the function $\mathcal{L}_{\lambda}\left(\cdot;\bm{y}\right)$. In the sequel, we denote by $\bm{L}_{\mathcal{G}}=\sum_{(i,j)\in\mathcal{E},i>j}\left(\bm{e}_{i}-\bm{e}_{j}\right)\left(\bm{e}_{i}-\bm{e}_{j}\right)^{\top}$ the (unnormalized) Laplacian matrix (Chung, 1997) associated with $\mathcal{G}$. For any matrix $\bm{A}$ we let

namely, the smallest eigenvalue when restricted to vectors orthogonal to $\bm{1}$.

Suppose that $p>\frac{c_{0}\log n}{n}$ for some sufficiently large constant $c_{0}>0$. Then on the event $\mathcal{A}_{0}$ as defined in (34), one has

Note that $\frac{e^{\theta_{i}}e^{\theta_{j}}}{(e^{\theta_{i}}+e^{\theta_{j}})^{2}}\leq\frac{1}{4}$. It follows immediately from the Hessian in (38) that

where $d_{\max}$ is the maximum vertex degree in the graph $\mathcal{G}$. In addition, on the event $\mathcal{A}_{0}$ we have $d_{\max}\leq 2np$, which completes the proof. ∎

Let $\mathcal{G}\sim\mathcal{G}_{n,p}$, and suppose that $p>\frac{c_{0}\log n}{n}$ for some sufficiently large constant $c_{0}>0$. Then one has

Note that $\lambda_{\min,\perp}\left(\bm{L}_{\mathcal{G}}\right)$ is exactly the spectral gap of the Laplacian matrix. See (Tropp, 2015, Sec 5.3.3) for the derivation of this lemma. ∎

By combining Lemma 9 with Lemma 10, we reach the following result.

Under the assumptions of Lemma 10, with probability exceeding $1-O\left(n^{-10}\right)$ one has

simultaneously for all $\bm{\theta}$ obeying $\left\|\bm{\theta}-\bm{\theta}^{*}\right\|_{\infty}\leq C$ for some $C\geq 0$.

2 Proof outline of Theorem 6

This subsection outlines the main steps for establishing Theorem 6.

Rather than directly resorting to the optimality condition, we adopt an algorithmic perspective to analyze the regularized MLE $\bm{\theta}$. Specifically, we consider the standard gradient descent algorithm that is expected to converge to the minimizer $\bm{\theta}$, and analyze the trajectory of this iterative algorithm instead. The algorithm is stated in Algorithm 2.

Notably, this gradient descent algorithm is not practical since the initial point is set to be $\bm{\theta}^{*}$. Nevertheless, it is helpful for analyzing the statistical accuracy of the regularized MLE $\bm{\theta}$. In what follows, we shall adopt a time-invariant step size rule:

Our proof can be divided into three steps:

establish — via standard optimization theory — that the output $\bm{\theta}^{T}$ of Algorithm 2 is sufficiently close to the regularized MLE $\bm{\theta}$, namely,

for $T=n^{5}$, where $C_{0}>0$ is some absolute constant;

use the leave-one-out argument to demonstrate that: the output $\bm{\theta}^{T}$ is close to the truth $\bm{\theta}^{*}$ in an entrywise fashion, i.e.

for some universal constant $C_{4}>0$. Combining this with (43) yields

the final step is to translate the perturbation bound on $\|\bm{\theta}-\bm{\theta}^{*}\|_{\infty}$ to $\|e^{\bm{\theta}}-e^{\bm{\theta}^{*}}\|_{\infty}$ as claimed in the theorem.

Before continuing, we single out an important fact that will be used throughout the proof.

Suppose $\bm{1}^{\top}\bm{\theta}^{*}=0$. Then we have $\bm{1}^{\top}\bm{\theta}^{t}=0$ for all $t\geq 0$.

3 Step I

The first step relies heavily on optimization theory, namely the theory of gradient descent on strongly convex and smooth functions.

It is seen that the sequence $\left\{\bm{\theta}^{t}\right\}_{t=1}^{\infty}$ converges geometrically fast to the regularized MLE $\bm{\theta}$, a property that is standard in convex optimization literature. This claim is summarized in the following lemma.

On the event $\mathcal{A}_{0}$ as defined in (34), one has

where $\rho=1-\frac{\lambda}{\lambda+np}$.

This result directly follows from the smoothness property (see Lemma 8), the trivial strong convexity of $\mathcal{L}_{\lambda}\left(\bm{\theta};\bm{y}\right)$ ($\nabla^{2}\mathcal{L}_{\lambda}\left(\bm{\theta};\bm{y}\right)\succeq\lambda\bm{I}_{n},\;\forall\bm{\theta}$), as well as the convergence property of the gradient descent algorithm (e.g. (Bubeck, 2015, Theorem 3.10)). ∎

A direct consequence of this convergence result and Fact 1 is that $\bm{1}^{\top}\bm{\theta}=0$ for the regularized MLE $\bm{\bm{\theta}}$.

We then control $\|\bm{\theta}^{0}-\bm{\theta}\|_{2}$. Recall that $\bm{\theta}^{0}=\bm{\theta}^{*}$, and we have:

On the event $\mathcal{A}_{2}$ as defined in (39), there exists some constant $c_{2}>0$ such that

The previous two claims taken together lead us to conclude that

for some constants $c_{3},c_{4},C_{0}>0$, $L\lesssim n^{5}$, and sufficiently large $T$ (recall that $T=n^{5}$). The above bounds are somewhat loose, but they suffice for our purpose. We then naturally obtain

as claimed. This finishes the first step of the proof.

4 Step II

where $\bm{\theta}^{0,\left(m\right)}=\bm{\theta}^{0}=\bm{\theta}^{*}$ and

Here, the leave-one-out loss function $\mathcal{L}_{\lambda}^{(m)}\left(\bm{\theta};\bm{y}\right)$ replaces all log-likelihood components involving the $m$-th item with their expected values (unconditional on $\mathcal{G}$). For any $1\leq m\leq n$, the auxiliary sequence $\left\{\bm{\theta}^{t,\left(m\right)}\right\}$ serves as a reasonably good proxy for $\left\{\bm{\theta}^{t}\right\}$, while remaining statistically independent of $\left\{y_{i,m}\mid(i,m)\in\mathcal{E}\right\}$.

Our proof in this step is inductive in nature. For the sake of clarity, we first list all induction hypotheses needed in our analysis:

where $C_{1},\cdots,C_{4}>0$ are some absolute constants. We aim to show that if the iterates at the $t$-th iteration — i.e. $\bm{\theta}^{t}$ and $\left\{\bm{\theta}^{t,\left(m\right)}\right\}_{1\leq m\leq n}$ — satisfy the induction hypotheses (46), then the $\left(t+1\right)$-th iterates continue to satisfy these hypotheses. Clearly, it suffices to justify (46) for all $0\leq t\leq T=n^{5}$.

Before we dive into the inductive arguments, there are a few direct consequences of (46) that are worth listing. We gather them in the next lemma.

Suppose the induction hypotheses (46) hold true for the $t$-th iteration, then there exist some universal constants $C_{5},C_{6}>0$ such that the following two bounds hold:

Note that the base case (i.e. the case for $t=0$) is trivially true due to the same initial points, namely, $\bm{\theta}^{0,\left(m\right)}=\bm{\theta}^{0}=\bm{\theta}^{*}$ for all $1\leq m\leq n$. We start with the first induction hypothesis (46a), which is supplied below.

Suppose the induction hypotheses (46) hold true for the $t$-th iteration, then with probability at least $1-O\left(n^{-10}\right)$, one has

as long as the step size obeys $0<\eta\leq\frac{1}{\lambda+np}$ and $C_{1}>0$ is sufficiently large.

The remaining induction steps are provided in the following lemmas.

Suppose the induction hypotheses (46) hold true for the $t$th iteration, then with probability at least $1-O\left(n^{-10}\right)$, one has

with the proviso that $0<\eta\leq\frac{1}{\lambda+np}$ and $C_{2}\gtrsim C_{6}+c_{\lambda}$.

Suppose the induction hypotheses (46) hold true for the $t$-th iteration, then with probability at least $1-O\left(n^{-10}\right)$, one has

as long as the step size obeys $0<\eta\leq\frac{1}{\lambda+np}$ and $C_{3}>0$ is sufficiently large.

Suppose the induction hypotheses (46) hold true for the $t$-th iteration, then with probability at least $1-O\left(n^{-10}\right)$, one has

Taking the union bound over $T=n^{5}$ iterations yields that with probability at least $1-O\left(n^{-5}\right)$,

which together with the conclusion in Step I results in

5 Step III

Toward this end, we observe that for each $1\leq m\leq n$,

as long as $\kappa^{2}\sqrt{\frac{\log n}{npL}}$ is small enough. This completes the proof of Theorem 6.

Acknowledgements

Y. Chen is supported in part by the grant ARO W911NF-18-1-0303 and by the Princeton SEAS innovation award. J. Fan is supported in part by NSF grants DMS-1662139 and DMS-1712591 and NIH grant 2R01-GM072611-13.

Appendix A Proof of Theorem 7

for some fixed constant $0<\epsilon<1$, then

Additionally, the inequality (54) comes from Pinsker’s inequality (Tsybakov, 2009, Lemma 2.5).

We then look at each term of (54) separately. To begin with, repeating the analysis in (Chen and Suh, 2015, Appendix B), we can demonstrate (using the independence assumption) that

where $\mathsf{Bern}(p)$ denotes the Bernoulli distribution with mean $p$. Upper bounding the KL divergence via $\chi^{2}$ divergence (see (Tsybakov, 2009, Lemma 2.7)), namely,

Put together the preceding bounds to reach

Appendix B Examples for the general κ𝜅\kappa setting

Recall that in Section 3, we introduce a new metric $\Delta_{K}^{*}$. The following three examples shed some light on the potential effectiveness of $\Delta_{K}^{*}$ as a fundamental information measure.

Case 1: $\kappa=O\left(1\right)$. Under this circumstance, it is easy to verify that

and hence all our preceding results for spectral method and regularized MLE for $\kappa=O(1)$ continue to hold with $\Delta_{K}$ replaced by $\Delta_{K}^{*}$. In this case, the new lower bound for sample complexity is slightly worse than the previous one (Theorem 2 in the main text) by a factor of $\log n$.

Case 2: Suppose there are 100 items with $w_{1}^{*}=\cdots=w_{5}^{*}=10$, $w_{6}^{*}=\cdots=w_{99}=5$ and $w_{100}^{*}=10^{-6}$. Our goal is to find the top-5 ranked items. Intuitively, the presence of the 100th item should not affect the hardness of top-5 ranking by much. This intuition is well captured by our new metric in (33) in the main text. Observe that

Since $w_{100}^{*}/w_{5}^{*}$ is exceedingly small ($10^{-7}$ in this example), it is easily seen that $\frac{w_{100}^{*}/w_{5}^{*}}{(1+w_{100}^{*}/w_{5}^{*})^{2}}\lesssim w_{100}^{*}/w_{5}^{*}$ is also extremely small, and hence $\left(\Delta_{K}^{*}\right)^{2}$ is not changed by much compared with the case when the 100th item is absent ($\left(\Delta_{K}^{*}\right)^{2}\approx 0.1107$ (resp. 0.1118) for the case when the $100$th item is present (resp. absent)). Similarly, consider the case when $w_{1}^{*}=10^{6}$, $w_{2}^{*}=\cdots=w_{5}^{*}=1$ and $w_{6}^{*}=\cdots=w_{100}^{*}=0.5$. As one can see, adding the first item will have little influence upon $\Delta_{K}^{*}$.

Case 3: Consider finding the top-5 items out of 100 items with $w_{1}^{*}=\cdots=w_{5}^{*}=10$, $w_{6}^{*}=\cdots=w_{10}^{*}=5$ and $w_{11}^{*}=\cdots w_{100}^{*}=10^{-6}$. The sample complexity needed for exact top-5 recovery will surely increase since the comparisons between $\left\{1,\cdots,10\right\}$ and $\left\{11,\cdots 100\right\}$ are, with high probability, not useful in determining the relative strength within the group $\left\{1,\cdots,10\right\}$. This is also reflected in the generalized separation measure $\Delta_{K}^{*}$. Recall that we have

For each $11\leq i\leq 100$, $\frac{w_{i}^{*}/w_{5}^{*}}{\left(1+w_{i}^{*}/w_{5}^{*}\right)^{2}}$ is exceedingly small. This makes $\left(\Delta_{K}^{*}\right)^{2}$ much smaller compared with the case when only $\left\{1,\cdots 10\right\}$ are present ($\left(\Delta_{K}^{*}\right)^{2}\approx 0.0118$ (resp. 0.1181) for the case when items $\{11,\cdots,100\}$ are present (resp. absent)). As a result, the required sample size increases accordingly.

Appendix C Proofs in Section 5

This section collects proofs of the theorems and lemmas that appear in Section 5.

Before moving on, we note that by Lemma 1 in the main text, the event

happens with probability at least $1-O(n^{-10})$. Throughout this section, we shall assume that we are on this event without explicitly referring to it each time. An immediate consequence is that $d_{\max}\leq d$ on this event.

The last term of the above identity can be further decomposed as

where we have used the fact that $\left(\bm{\pi}-\hat{\bm{\pi}}\right)^{\top}\bm{1}\bm{\pi}^{*\top}=\bm{0}$. Combining (55) and (56) we get

which together with a little algebra gives

C.2 Proof of Theorem 9

where (i) is a consequence of Lemma 5, (ii) follows from the relationship between $\|\cdot\|_{2}$ and $\|\cdot\|_{\bm{\pi}^{*}}$, and (iii) follows as long as one can justify that

Therefore, the rest of the proof is devoted to establishing (57). To simplify the notations hereafter, we denote $\bm{\Delta}:=\bm{P}-\bm{P}^{*}$. In fact, it is easy to check that for any $i\neq j$,

In view of Hoeffding’s inequality (Lemma 18), one has, when conditional on $\mathcal{G}$, that

for some constant $c>0$. By choosing $t\asymp\sigma^{2}\sqrt{n\log n}$, we see that with probability at least $1-O(n^{-10})$,

The same upper bounds can be derived for other terms using the same arguments. We have thus established (57) by recognizing that $d\gtrsim np$.

C.3 Proof of Lemma 2

where $\left(\text{i}\right)$ follows from the fact that $\bm{P}$ and $\bm{P}^{*}$ are both probability transition matrices. By Lemma 18, one can derive

When $d_{\max}\leq d$, the right hand side is bounded by $2\exp\left(\frac{-Ldt^{2}}{2\|\bm{\pi}^{*}\|_{\infty}^{2}}\right)$. Hence

C.4 Proof of Lemma 3

Similar to the proof of Lemma 2 above, by applying Lemma 18 to the quantity

On the other hand, when the event $\mathcal{A}_{0}$ (defined in Lemma 1 in the main text) happens, we have $d_{\min}\geq np/2$ and for all $1\leq m\leq n$,

Combining these two pieces completes the proof.

C.5 Proof of Lemma 4

First, by the relationship between $\|\cdot\|_{2}$ and $\|\cdot\|_{\bm{\pi}^{*}}$, we have

where $\pi^{*}_{\min}:=\min_{i}\pi_{i}^{*}$. Invoking Theorem 8, we obtain

where we define $\gamma:=1-\max\left\{\lambda_{2}(\bm{P}^{*}),-\lambda_{n}\left(\bm{P}^{*}\right)\right\}-\|\bm{P}-\bm{P}^{*}\|_{\bm{\pi}^{*}}$ and $\pi^{*}_{\max}:=\max_{i}\pi_{i}^{*}$.

To facilitate the analysis of $\|\bm{\pi}^{\left(m\right)\top}(\bm{P}^{(m)}-\bm{P})\|_{2}$, we introduce another Markov chain with transition probability matrices $\bm{P}^{(m),\mathcal{G}}$, which is also a leave-one-out version of the transition matrix $\bm{P}$. Similar to $\bm{P}^{(m)}$, $\bm{P}^{(m),\mathcal{G}}$ replaces all the transition probabilities involving the $m$-th item with their expected values (conditional on $\mathcal{G}$). Concretely, for $i\neq j$,

to make $\bm{P}^{(m),\mathcal{G}}$ a valid probability transition matrix. Hence by the triangle inequality, we see that

The next step is then to bound $J_{1}^{m}$ and $J_{2}^{m}$ separately.

For $J_{1}^{m}$, similar to (60), one has

where $\left(\text{i}\right)$ comes from the fact that $\bm{P}_{\cdot m}^{(m),\mathcal{G}}=\bm{P}_{\cdot m}^{*}$. Recognizing that $\bm{\pi}^{\left(m\right)}$ is statistically independent of $\left\{y_{j,m}\right\}_{j\neq m}$, by Hoeffding’s inequality in Lemma 18, we get

In addition by Hoeffding’s inequality in Lemma 18, we have

with probability at least $1-O(n^{-5})$. As a consequence,

where $(\text{i})$ comes from the fact that $d_{\max}\leq d$.

Regarding $J_{2}^{m}$, we invoke the identity $\bm{\pi}^{*\top}(\bm{P}^{(m)}-\bm{P}^{(m),\mathcal{G}})=\bm{0}$ to get

Recognizing that \big{|}P_{j,m}^{(m)}-P_{j,m}^{(m),\mathcal{G}}\big{|}\leq\frac{2}{d} for $\left(j,m\right)\in\mathcal{E}$ and \big{|}P_{j,m}^{(m)}-P_{j,m}^{(m),\mathcal{G}}\big{|}\leq\frac{p}{d} for $\left(j,m\right)\notin\mathcal{E}$, we have

Given that $P_{m,j}^{(m)}-P_{m,j}^{(m),\mathcal{G}}=\frac{y_{m,j}^{*}}{d}(p-\operatorname{\mathds{1}}_{(m,j)\in\mathcal{E}})$, we have

Since $\|\bm{\xi}^{(m)}\|_{\infty}\leq\frac{1}{d}\|\bm{\pi}^{(m)}-\bm{\pi}^{*}\|_{\infty}$ and $\|\bm{\xi}^{(m)}\|_{2}\leq\frac{1}{d}\|\bm{\pi}^{(m)}-\bm{\pi}^{*}\|_{2}$, Lemma 19 implies that

with high probability. The same bound holds for $J_{4}^{m}$. Combine (63), (65) and (66) to arrive at

where $(\text{i})$ holds as long as $np\gamma^{2}\geq c\kappa\log n$ for $c$ sufficiently large. The triangle inequality

C.6 Proof of Lemma 6

We can further decompose $I_{4}^{m}$ into

where $\left(\text{i}\right)$ comes from the Cauchy-Schwarz inequality, $\left(\text{ii}\right)$ follows from the choice $d=c_{d}np\geq 2np$ and $\left(\text{iii}\right)$ results from the triangle inequality. By Theorem 9, with high probability we have

When it comes to the fluctuation term, one can write

Since $\|\bm{\beta}^{(m)}\|_{2}\leq\frac{1}{d}\|\bm{\pi}^{(m)}-\bm{\pi}^{*}\|_{2}$ and $\|\bm{\beta}^{(m)}\|_{\infty}\leq\frac{1}{d}\|\bm{\pi}^{(m)}-\bm{\pi}^{*}\|_{\infty}$, one can apply Lemma 19 to derive

with high enough probability. The bounds (68) and (69) taken together complete the proof.

Appendix D Proofs in Section 6

This section gathers the proofs of the lemmas in Section 6.

This implies that with high probability (note that the randomness comes from $\mathcal{G}$),

where the last relation holds because of the facts that $\left\|\bm{\theta}^{*}\right\|_{2}\leq\sqrt{n}\log\kappa$, $\lambda\asymp\frac{1}{\log\kappa}\sqrt{\frac{np\log n}{L}}$, and

D.2 Proof of Lemma 9

where the last relation holds since $\left(1+e^{-x}\right)^{2}\leq 4$ for all $x\geq 0$. From our assumption, we see that for all $1\leq i,j\leq n$,

which relies on the fact that $\theta_{\max}^{*}-\theta_{\min}^{*}\leq\log\kappa$. This allows one to justify that

D.3 Proof of Fact 1

By $\bm{\theta}^{0}=\bm{\theta}^{*}$ and $\bm{1}^{\top}\bm{\theta}^{*}=0$, the statement trivially holds true for $t=0$. Suppose it is true for some $t\geq 0$. Then

where the equalities (i) and (ii) follow from the fact that $\bm{1}^{\top}\bm{\theta}^{t}=0$, whereas the last identity (iii) arises from the gradient expression (37) and the simple fact that $\bm{1}^{\top}(\bm{e}_{i}-\bm{e}_{j})=0$ for any $i$ and $j$. This completes the whole proof.

D.4 Proof of Lemma 12

It follows from the optimality of $\bm{\theta}$ as well as the mean value theorem that

On the event $\mathcal{A}_{2}=\left\{\left\|\nabla\mathcal{L}_{\lambda}\left(\bm{\theta}^{*};\bm{y}\right)\right\|_{2}\lesssim\sqrt{\frac{n^{2}p\log n}{L}}\right\}$ and in the presence of the choice $\lambda\asymp\frac{1}{\log\kappa}\sqrt{\frac{np\log n}{L}}$, we obtain $\left\|\bm{\theta}-\bm{\theta}^{*}\right\|_{2}\leq c_{2}\log\kappa\sqrt{n}$ for some constant $c_{2}>0$.

D.5 Proof of Lemma 13

In regard to the first consequence, one can apply the triangle inequality to show

as long as $C_{5}\geq C_{4}+C_{3}$. Similarly, for the second one, we have

D.6 Proof of Lemma 14

In view of the gradient update rule (41), we have

where we denote $\bm{\theta}\left(\tau\right):=\bm{\theta}^{*}+\tau\left(\bm{\theta}^{t}-\bm{\theta}^{*}\right)$. Here, the last identity results from the fundamental theorem of calculus (Lang, 1993, Chapter XIII, Theorem 4.2). Let $\theta_{\max}\left(\tau\right):=\max_{i}\theta_{i}(\tau)$ and $\theta_{\min}\left(\tau\right):=\min_{i}\theta_{i}(\tau)$. Combining the induction hypothesis (46d) with the definition of $\bm{\theta}\left(\tau\right)$, one can see that for all $0\leq\tau\leq 1$,

for any sufficiently small $\epsilon>0$, as long as

This together with Lemma 8 and Corollary 1 reveals that for any $0\leq\tau\leq 1$,

Since $\mathbf{1}^{\top}(\bm{\theta}^{t}-\bm{\theta}^{*})=0$ , the first term on the right hand side of (74) is controlled by

Substitute (73) into the above inequality to reach

Substitute (75) back to (74) and use the induction hypothesis (46d) to conclude that

for some constants $C,C_{1}>0$. Here, the second inequality makes use of the facts that $\left\|\nabla\mathcal{L}_{\lambda}\left(\bm{\theta}^{*}\right)\right\|_{2}\lesssim\sqrt{\frac{n^{2}p\log n}{L}}$ (see Lemma 7). The last line holds with the proviso that $C_{1}>0$ is sufficiently large.

D.7 Proof of Lemma 15

Consider any $m$ ($1\leq m\leq n$). According to the gradient update rule (44), one has

where the last line follows by the construction of $\mathcal{L}_{\lambda}^{\left(m\right)}$. Apply the mean value theorem to obtain