Coordinate Descent with Arbitrary Sampling II: Expected Separable Overapproximation

Introduction

Coordinate descent methods have been popular with practitioners for many decades due to their inherent conceptual simplicity and ease with which one can produce a working code. However, up to a few exceptions , they have been largely ignored in the optimization community until recently when a renewed interest in coordinate descent was sparked by several reports of their remarkable success in certain applications . Additional and perhaps more significant reason behind the recent flurry of research activity in the area of coordinate descent comes from breakthroughs in our theoretical understanding of these methods through the introduction of randomization in the iterative process . Traditional variants of coordinate descent rely on cyclic or greedy rules for the selection of the next coordinate to be updated.

It has recently become increasingly clear that the design and complexity analysis of randomized coordinate descent methods is intimately linked with and can be better understood through the notion of expected separable overapproximation (ESO) and . This refers to an inequality involving the objective function and the sampling (a random set valued mapping describing the law with which subsets of coordinates are selected at each iteration), capturing in a compact way certain smoothness properties of the function in a random subspace spanned by the sampled coordinates.

A (coordinate) sampling $\hat{S}$ is a random set-valued mapping with values being subsets of $[n]\stackrel{{\scriptstyle\text{def}}}{{=}}\{1,2,\dots,n\}$. It will be useful to write

We will compactly write $(f,\hat{S})\sim ESO(v)$.

Instead of the above general definition, it will be useful to the reader to instead think about the form of this inequality in the simple case when $f(x)=\|Ax\|^{2}$, where $\|\cdot\|$ is the L2 norm, and $x=0$. Letting $A=[A_{1},\dots,A_{n}]$, in this case inequality (2) takes the form

The ESO inequality is of key importance for randomized coordinate descent methods for several reasons:

The parameters $v=(v_{1},\dots,v_{n})$ for which ESO holds are neededAll existing parallel coordinate coordinate descent methods for which a complexity analysis has been performed are designed with fixed stepsizes. Designing a line-search procedure in such a setup is a nontrivial task, and to the best of our knowledge, only a single paper in the literature deals with this issue . Certainly, properly designed line search has the potential to improve the practical performance of these methods. to run coordinate descent. Indeed, they are used to set the stepsizes to a suitable value.

The size of these parameters directly influences the complexity of the method (see Table 1).

There are problems for which updating more coordinates in each iteration, as opposed to updating just one, may not lead to fewer iterations (which suggests that perhaps the resources should be instead utilized in some other way). Whether this happens or not can be understood through a careful study of the complexity result and its dependence, through the vectors $p$ and $v$, on the number of coordinates updated in each iteration .

The ESO assumption is generic in the sense that as soon as function $f$ and sampling $\hat{S}$ satisfy it, the complexity result follows. This leads to a natural dichotomy in the study of coordinate descent: i) the search for new variants of coordinate descent (e.g., parallel, accelerated, distributed) and study of their complexity under the ESO assumption, and ii) the search for pairs $(f,\hat{S})$ for which one can compute $v$ such that $(f,\hat{S})\sim ESO(v)$. Our current study follows this dichotomy: in we deal with the algorithmic and complexity aspects, and in this paper we deal with the ESO aspect.

2 Complexity of coordinate descent

For instance, the NSync boundComplexity of NSync depends on the initial ($x_{0}$) and optimal ($x_{*}$) points, but have hidden this dependence. The full bound is obtained by replacing $\log(1/\epsilon)$ by $\log((f(x_{0})-f(x_{*}))/\epsilon)$, where $f$ is the objective function. in Table 1 applies to the problem of unconstrained minimization of a smooth strongly convex function. It was in where the general form of the ESO inequality used in this paper was first mentioned and used to derive a complexity result for a coordinate descent method with arbitrary sampling.

The Quartz algorithm , on the other hand, applies to a much more serious problem – a problem of key importance in machine learning. In particular, it applies to the regularized empirical risk minimization problem, where the loss functions are convex and have Lipschitz gradients and the regularizer is strongly convex and possibly nonsmooth. Coordinate ascent is applied to the dual of this problem, and the bound appearing in Table 1 applies to the duality gapComplexity of Quartz depends on an initial pair of primal and dual vectors; we have omitted this dependence from the table. The full complexity result is obtained by replacing $\log(1/\epsilon)$ by $\log(\Delta_{0}/\epsilon)$, where $\Delta_{0}$ is the difference between the primal and dual function values for a pair of (primal and dual) starting points..

The APPROX method was first proposed in and then generalized to an arbitrary sampling (among other things) in . In its accelerated variant it enjoys a $O(1/\sqrt{\epsilon})$ rate, whereas it’s non-accelerated variant has a slower $O(1/\epsilon)$ rate. Again, the complexity of the method explicitly depends on the vector of probabilities $p$ and the ESO parameter $v$.

3 Historical remarks

4 Contributions

We now briefly list the contributions of this work.

ESO inequalities were previously established for special classes of samplings only, almost invariably for uniform samplings , and often using seemingly disparate approaches. We give the first systematic study of ESO inequalities for arbitrary samplings.

We recover existing ESO results by applying our general technique.

Our approach to deriving ESO inequalities is via the study of random principal submatrices of a positive semidefinite matrix. In particular, we give bounds on the largest eigenvalue of the mean of the random submatrix. This may be of independent interest.

5 Outline of the paper

Our paper is organized as follows. In Section 2 we describe the class of functions ($f$) we consider in this paper and briefly establish some basic terminology related to samplings ($\hat{S}$). In Section 3 we study probability matrices associated with samplings ($\mathbf{P}(\hat{S})$), in Section 4 we study eigenvalues of these probability matrices ($\lambda(\mathbf{P}(\hat{S}))$ and $\lambda^{\prime}(\mathbf{P}(\hat{S}))$) and in Section 5 we design a general technique for computing parameter $v=(v_{1},\dots,v_{n})$ for which the ESO inequality holds (i.e., for which $(f,\hat{S})\sim ESO(v)$). We illustrate the use of these techniques in Section 5.4 and conclude with Section 7.

Functions and samplings

Recall that in the paper we are concerned with establishing inequality (2) which we succinctly write as $(f,\hat{S})\sim ESO(v)$. In Section 2.1 we describe the class of functions $f$ we consider in this paper and in Section 2.2 we briefly review several elementary facts related to samplings.

In the subsequent text, we shall often refer to the set of columns of $\mathbf{A}$ for which the entry in the $j$-th row of $\mathbf{A}$ is nonzero:

Assumption 2.1 holds for many functions of interest in optimization and machine learning. Coordinate descent methods for functions $f$ explicitly required to satisfy Assumption 2.1 were studied in .

The following simple observation will help us relate the above assumption with standing assumptions considered in various papers on randomized coordinate descent methods.

It remains to add these inequalities for $j=1,\dots,s$. ∎

By $\mathbf{I}$ we denote the $n$-by-$n$ identity matrix and for $S\subseteq[n]$ we will use the notation $\mathbf{I}_{[S]}$ for the $n$-by-$n$ matrix obtained from $\mathbf{I}$ by retaining elements $\mathbf{I}_{ii}$ for which $i\in S$ and zeroing out all other elements.

We now apply Proposition 2.1 to several special cases:

Partial separability. Let $d=n$ and $\mathbf{M}_{j}=\mathbf{I}_{[C_{j}]}$, where for each $j$, $C_{j}\subseteq[n]$. Then $f$ is of the form

That is, $\phi_{j}$ depends on coordinates of $x$ belonging to set $C_{j}$ only. By Proposition 2.1, $f$ satisfies (3), where $\mathbf{A}$ is the $n$-by-$n$ diagonal matrix given by

Functions of the form (6) (i.e., partially separable functions) were considered in the context of parallel coordinate descent methods in . However, in the authors only assume the sum $f$ to have a Lipschitz gradient (which is more general, but somewhat complicates the analysis), whereas we assume that all component functions $\{\phi_{j}\}_{j}$ have Lipschitz gradient.

Linear transformation of variables. Let $s=1$. Then $f$ is of the form

By Proposition 2.1, $f$ satisfies (3), where $\mathbf{A}$ is given by

A functions of the form (7) appears in the dual problem of the standard primal-dual formulation to which stochastic dual coordinate ascent methods are applied [27, MinibatchASDCA, 33, 10, 18].

By Proposition 2.1, $f$ satisfies (3), with $\mathbf{A}$ given by

Functions of the form (8) play an important role in the design of efficiently implementable accelerated coordinate descent methods . These functions also appear in the primal problem of the standard primal-dual formulation to which stochastic dual coordinate ascent methods are applied.

2 Samplings

As defined in the introduction, by sampling we mean a random set-valued mapping with values in $2^{[n]}$ (the set of subsets of $[n]$).

Of key importance in this paper are elementary samplings, defined next.

Let ${\cal P}_{1},\dots,{\cal P}_{c}$ be a partition of $\{1,2,\dots,n\}$ such that $|{\cal P}_{l}|=s$ for all $l$. That is, $sc=n$. Now let $\hat{S}_{1},\dots,\hat{S}_{c}$ be independent $\tau$-nice samplings from ${\cal P}_{1},\dots,{\cal P}_{c}$, respectively. Then the sampling

is called $(c,\tau)$-distributed sampling.

The $\tau$-nice sampling arises as a special case of the $(c,\tau)$-distributed sampling (for $c=1$) which we define next.

Sampling $\hat{S}$ is called $\tau$-nice if it picks only subsets of $[n]$ of cardinality $\tau$, uniformly at random. More formally, it is defined by

Operations with samplings.

We now define several basic operations with samplings (convex combination, intersection and restriction).

Let $\hat{S}_{1},\dots,\hat{S}_{k}$ be samplings and let $q_{1},\dots,q_{k}$ be nonnegative scalars summing to 1. By $\sum_{t=1}^{k}q_{t}\hat{S}_{t}$ we denote the sampling obtained as follows: we first pick $t\in\{1,\dots,k\}$, with probability $q_{t}$, and then sample according to $\hat{S}_{t}$. More formally, $\hat{S}$ is defined as follows:

Note that (11) indeed defines a sampling, since

Each sampling is a convex combination of elementary samplings. Indeed, for each $\hat{S}$ we have

We now show that each doubly uniform sampling arises as a convex combination of $\tau$-nice samplings.

Let $\hat{S}$ be a doubly uniform sampling and let $\hat{S}_{\tau}$ be the $\tau$-nice sampling, for $\tau=0,1,\dots,n$. Then

where the last equality follows from the definition of doubly uniform and $\tau$-nice samplings. The statement then follows from (11) (i.e., by definition of convex combination of samplings). ∎

It will be useful to define two more operations with samplings; intersection and restriction.

For two samplings $\hat{S}_{1}$ and $\hat{S}_{2}$ we define the intersection $\hat{S}\stackrel{{\scriptstyle\text{def}}}{{=}}\hat{S}_{1}\cap\hat{S}_{2}$ as the sampling for which:

Let $\hat{S}$ be a sampling and $J\subseteq[n]$. By restriction of $\hat{S}$ to $J$ we mean the sampling $\hat{E}_{J}\cap\hat{S}$. By abuse of notation we will also write this sampling as $J\cap\hat{S}$.

Graph sampling.

Let $G=(V,E)$ be an undirected graph with $|V|=n$ vertex and $(i,i^{\prime})$ be an edge in $E$ if and only if there is $j\in[m]$ such that $\{i,i^{\prime}\}\subseteq J_{j}$. If $S$ is an independent set of graph $G$, then necessarily

Denote by $\mathcal{T}$ the collection of all independent sets of the graph $G$. We now define the graph sampling as follows:

Let $\hat{S}$ be a graph sampling. In view of (12), for some nonnegative constants $q_{S}$ adding up to 1:

Let $X_{1},\dots,X_{\tau}$ be a partition of $[n]$, i.e.,

The product sampling $\hat{S}$ is obtained by choosing $S\in{\cal S}$, uniformly at random; that is, via:

A similar sampling was first considered in [18, Section 3.3] with an additional group separability assumption on the partition $X_{1},\dots,X_{\tau}$, which can be equivalently stated as:

In other words, it is both a product sampling and graph sampling. Note that in Definition 2.8 we do not make any assumption on the partition. Also, the product sampling is a nonuniform sampling as long as all the sets $X_{l}$ do not have the same cardinality, which occurs necessarily if $\tau$, representing the number of processors, is not divisible by $n$.

Probability matrix associated with a sampling

In this section we define the notion of a probability matrix associated with a sampling. As we shall see in later sections, this matrix encodes all information about $\hat{S}$ which is relevant for development of ESO inequality.

With each sampling $\hat{S}$ we associate an $n$-by-$n$ “probability matrix” $\mathbf{P}=\mathbf{P}(\hat{S})$ defined by

We shall write $\mathbf{P}(\hat{S})$ when it is important to indicate which sampling is behind the probability matrix, otherwise we simply write $\mathbf{P}$.

Using the notation we have just established, probability matrices of elementary samplings are given by

We now establish a simple but particularly insightful result, leading to many useful identities.

The set of probability matrices is the convex hull of the probability matrices corresponding to elementary samplings.

$\mathbf{P}(\hat{S})\succcurlyeq 0$ for each $\hat{S}$.

Since multiplying a matrix in the Hadamard sense by a fixed matrix is a linear operation,

Identity (20) follows from (19) by setting $\mathbf{M}=\mathbf{E}$:

Finally, (22) (resp. (23)) follows from (20) (resp. (21)) by setting $h=e$. ∎

2 Operations with samplings

We now give formulae for the probability matrix of the sampling arising as a convex combination, intersection or a restriction, in terms of the probability matrices of the constituent samplings.

We have seen in (12) that each sampling is a convex combination of elementary samplings. In view of Theorem 3.1, the probability matrices of the samplings are related the same way:

More generally, as formalized in the following lemma, the probability matrix of a convex combination of samplings is equal to the convex combination of the probability matrices of these samplings.

Let $\hat{S}_{1},\dots,\hat{S}_{k}$ be samplings and $q_{1},\dots,q_{k}$ be non-negative scalars summing up to 1. Then

Let $\hat{S}$ be the convex combination of samplings $\hat{S}_{1},\dots,\hat{S}_{k}$ and fix any $i,j\in[n]$. By definition,

Intersection of samplings.

The probability matrix of the intersection of two independent samplings is equal to the Hadamard product of the probability matrices of these samplings. This is formalized in the following lemma.

Let $\hat{S}_{1},\hat{S}_{2}$ be independent samplings. Then

Restriction.

By Lemma 3.2, the probability matrix of the restriction of arbitrary sampling $\hat{S}$ to $J\subseteq[n]$ is given by (we give several alternative ways of writing the result):

Note that $\mathbf{P}(J\cap\hat{S})$ is the matrix obtained from $\mathbf{P}(\hat{S})$ by keeping only elements $i,j\in J$ and zeroing out all the rest. Furthermore, by combining the formulae derived above, we get

3 Probability matrix of special samplings

The probability matrix of the $(c,\tau)$-distributed samplings is computed in the following lemma.

Let $\hat{S}$ be the ($c,\tau$)-distributed sampling associated with the partition $\{\mathcal{P}_{1},\dots,\mathcal{P}_{c}\}$ of $[n]$ such that $s=|\mathcal{P}_{l}|$ for $l\in[c]$ (see Definition 2.2). Then

Note that $\mathbf{B}$ is the 0-1 matrix with $\mathbf{B}_{ij}=1$ if and only if $i,j$ belong to the same partition.

Let $\mathbf{P}=\mathbf{P}(\hat{S})$. It is easy to see that

As a corollary of the above in the $c=1$ case we obtain the probability matrix of the $\tau$-nice sampling:

Fix $1\leqslant\tau\leqslant n$ and let $\hat{S}$ be the $\tau$-nice sampling. Then

where $\beta=(\tau-1)/\max(n-1,1)$. If $\tau=0$, then $\mathbf{P}(\hat{S})$ is the zero matrix.

For $\tau\geqslant 1$ this follows from Lemma 3.3 in the special case when $c=1$ (note that ${\cal P}_{1}=[n],s=n$ and $\mathbf{B}=\mathbf{E}$). ∎

Finally, we compute the probability matrix of a doubly uniform sampling.

Largest eigenvalues of the probability matrix

For an $n\times n$ positive semidefinite matrix $\mathbf{M}$ we denote by $\lambda(\mathbf{M})$ the largest eigenvalue of $\mathbf{M}$:

Note that $1\leqslant\lambda^{\prime}(\mathbf{M})\leqslant n$.

In this section we study (standard and normalized) largest eigenvalue of the probability matrix associated with a sampling:

Recall that by Theorem 3.1, $\mathbf{P}(\hat{S})$ is positive semidefinite for each sampling $\hat{S}$. For convenience, we write $\lambda(\hat{S})$ (resp. $\lambda^{\prime}(\hat{S})$) instead of $\lambda(\mathbf{P}(\hat{S}))$ (resp. $\lambda^{\prime}(\mathbf{P}(\hat{S}))$). We study these quantities since, as we will show in later sections, they are useful in computing parameter $v=(v_{1},\dots,v_{n})$ for which ESO holds.

In the case of elementary samplings the situation is simple. Indeed, for any $J\subseteq[n]$, we have

Since $\mathbf{P}(\hat{E}_{J})=\mathbf{E}_{[J]}$ and $\operatorname{Diag}(\mathbf{E}_{[J]})=\mathbf{I}_{[J]}$, (39) can equivalently be written as

2 Bounds for arbitrary samplings

In the first result of this section we give sharp bounds for $\lambda^{\prime}(\hat{S})$ for arbitrary sampling $\hat{S}$.

Upper bound. If $\tau$ is a constant such that $|\hat{S}|\leqslant\tau$ with probability 1, then $\lambda^{\prime}(\hat{S})\leqslant\tau.$

Identity. If $|\hat{S}|=\tau$ with probability 1, then $\lambda^{\prime}(\hat{S})=\tau$.

where the last inequality holds since $\operatorname{Tr}(\mathbf{P})$ is upper bounded by the sum of all elements of $\mathbf{P}$. It remains to apply identities (22) and (23).

In view of (12), we can represent $\hat{S}$ as a convex combination of elementary samplings:

The result follows by combining the upper and lower bounds. ∎

In the next result we study the quantity $\lambda(\hat{S})$.

Lower and upper bounds. For any sampling $\hat{S}$ we have

Sharper upper bound. If $\hat{S}$ is uniform and $|\hat{S}|\leqslant\tau$ with probability one, then the upper bound can be improved to

Identity. If $\hat{S}$ is uniform and $|\hat{S}|=\tau$ with probability one, then

By combining (38) and Theorem 4.1 (ii) we obtain:

The result follows by combining the lower bound from (i) with the upper bound in (ii).∎

3 Bounds for restrictions of selected samplings

In this part we study the normalized eigenvalue associated with the restriction of a few selected samplings (or families of samplings). In particular, we first give a (necessarily rough) bound that holds for arbitrary samplings, followed by a bound for the $(c,\tau)$-distributed sampling and the $\tau$-nice sampling (both are specific uniform samplings). Finally, we give a bound for the family of doubly uniform samplings.

Let $\hat{S}$ be an arbitrary sampling and let $\tau$ be such that $|\hat{S}|\leqslant\tau$ with probability 1. Then for all $\emptyset\neq J\subseteq[n]$, we have

Note that $|J\cap\hat{S}|\leqslant\min\{|J|,\tau\}$ with probability 1. We only need to apply the upper bound in Theorem 4.1 to the restriction sampling $J\cap\hat{S}$. ∎

We now proceed to the $(c,\tau)$-distributed sampling (recall Definition 2.2).

Let $\hat{S}$ be the ($c,\tau$)-distributed sampling associated with a partition $\{\mathcal{P}_{1},\dots,\mathcal{P}_{c}\}$ of $[n]$ such that $s=|\mathcal{P}_{l}|$ for $l\in[c]$. Fix arbitrary $\emptyset\neq J\subseteq[n]$ and let $\omega^{\prime}$ be the number of sets $\mathcal{P}_{l}$ which have a nonempty intersection with $J$; that is, let $\omega^{\prime}\stackrel{{\scriptstyle\text{def}}}{{=}}|\{l:J\cap\mathcal{P}_{l}\neq\emptyset\}|$. Then

By applying Lemma 3.2 and Lemma 3.3, we get

where the inequality is an application of the Cauchy-Schwartz inequality. It follows that

Finally, note that $\operatorname{Diag}(\mathbf{P}(\hat{E}_{J}))\circ\operatorname{Diag}(\mathbf{P}(\hat{S}))=\operatorname{Diag}(\mathbf{P}(\hat{E}_{J})\circ\mathbf{P}(\hat{S}))\overset{\eqref{eq:98shsguy8658}}{=}\operatorname{Diag}(\mathbf{P}(J\cap\hat{S}))$. ∎

We now specialize the above result to the $c=1$ case, obtaining a formula for $\lambda^{\prime}(J\cap\hat{S})$ in the case when $\hat{S}$ is the $\tau$-nice sampling (recall Definition 2.3).

Let $\hat{S}$ be the $\tau$-nice sampling. Then for all $\emptyset\neq J\subseteq[n]$,

Let $\emptyset\neq J\subseteq[n]$. Since $\tau$-nice sampling is the $(1,\tau)$-distributed sampling, by applying Proposition 4.2 we get:

Next, by direct calculation we can verify that

which together with the lower bound established in Theorem 4.1 yields:

Note that (48) is much better (i.e., smaller) than the right hand side in (43). This is to be expected as the bound (43) applies to all samplings (which have size at most $\tau$ with probability 1).

Finally, we give a bound on the normalized largest eigenvalue of the restriction of a doubly uniform sampling.

Expected Separable Overapproximation

In this section we develop a general technique for computing parameters $v=(v_{1},\dots,v_{n})$ for which the ESO inequality (2) holds.

The reason for defining and studying probability matrices $\mathbf{P}(\hat{S})$ is motivated by the following result, which for functions satisfying Assumption 2.1 reduces the ESO Assumption $(f,\hat{S})\sim ESO(v)$ to the problem of bounding the Hadamard product of the probability matrix $\mathbf{P}(\hat{S})$ and the data matrix $\mathbf{A}^{\top}\mathbf{A}$ from above by a diagonal matrix. Note that because $\mathbf{P}(\hat{S})\succcurlyeq 0$, in view of (50),the Hadamard product $\mathbf{P}(\hat{S})\circ\mathbf{A}^{\top}\mathbf{A}$ is positive semidefinite.

Let us substitute $h\leftarrow h_{[\hat{S}]}$ into (3) and take expectation in $\hat{S}$ of both sides. Applying (19), we obtain:

We next focus on the problem of finding vector $v$ for which (51) holds. The following direct consequence of (50) will be helpful in this regard:

In particular, (53) can be used to establish the first part of the following useful lemma.

If $\mathbf{M}_{1}\succcurlyeq 0$ and $\mathbf{M}_{2}\succcurlyeq 0$, then

By definition, $\mathbf{M}_{2}\preccurlyeq\lambda^{\prime}(\mathbf{M}_{2})\operatorname{Diag}(\mathbf{M}_{2})$, which together with (53) implies:

Applying the same reasoning to the matrix $\mathbf{M}_{1}$ we obtain: $\mathbf{M}_{1}\circ\mathbf{M}_{2}\preccurlyeq\lambda^{\prime}(\mathbf{M}_{1})\operatorname{Diag}(\mathbf{M}_{1}\circ\mathbf{M}_{2}).$ Combining the two results, we obtain (54). Inequality (55) follows from:

2 ESO I: no coupling between the sampling and data

By applying Lemma 5.2, Eq (54), to $\mathbf{M}_{1}=\mathbf{P}(\hat{S})$ and $\mathbf{M}_{2}=\mathbf{A}^{\top}\mathbf{A}$, we obtain a formula for $v$ satisfying (51).

Let $f$ satisfy Assumption 2.1 and let $\hat{S}$ be an arbitrary sampling. Then $(f,\hat{S})\sim ESO(v)$ for $v=(v_{1},\dots,v_{n})$ defined by

Let $\mathbf{P}=\mathbf{P}(\hat{S})$. To establish the main statement, it is sufficient to apply Lemma 5.1 and Lemma 5.2 and note that for $v$ defined by (56), $\operatorname{Diag}(v\circ p)=\min(\lambda^{\prime}(\mathbf{P}),\lambda^{\prime}(\mathbf{A}^{\top}\mathbf{A}))\operatorname{Diag}(\mathbf{P}\circ\mathbf{A}^{\top}\mathbf{A})$. ∎

If for some $\tau$, $|\hat{S}|\leqslant\tau$ with probability 1, then in view of Theorem 4.1, we have $\lambda^{\prime}(\mathbf{P}(\hat{S}))\leqslant\tau$. Furthermore,

Hence, in view of Lemma 5.1, we can pick the ESO parameter conservatively as follows:

An ESO inequality with $v_{i}$ similar to (57) was established in , but for a different class of functions ($\omega$-partially separable functions: functions expressed as a sum of functions each of which depends on at most $\omega$ coordinates) and uniform samplings only. Indeed, the bound established therein for arbitrary uniform samplings uses $v_{i}=\min\{\tau,\omega\}L_{i}$, where $\omega$ is the degree of separability of $f$ and $L_{i}$ is the Lipschitz constant of $\nabla f$ associated with coordinate $i$. In our setting, $\omega=\max_{j}|J_{j}|$ and $L_{i}$ corresponds to $\sum_{j}\mathbf{A}_{ji}^{2}$. Hence, (57) could be seen as a generalization of the ESO bound in to arbitrary samplings.

Note that computation of the normalized eigenvalue $\lambda^{\prime}(\mathbf{A}^{\top}\mathbf{A})$ could be time-consuming, and would require a number of passes through the data prior to running a coordinate descent method, which may be prohibitive. In the next section we follow a different approach, one in which this issue is avoided. The main idea is to decompose $\mathbf{A}^{\top}\mathbf{A}$ as a sum of the rank one matrices $\mathbf{A}_{j:}^{\top}\mathbf{A}_{j:}$ and then bound each term $\mathbf{P}(\hat{S})\circ\mathbf{A}_{j:}^{\top}\mathbf{A}_{j:}$ separately.

3 ESO II: coupling the sampling with data

In this section we use a different strategy for satisfying (51). We first write

where $\mathbf{A}_{j:}$ denote the $j$th row vector of matrix $\mathbf{A}$ and then bound each term in the last sum individually. Recall the definition of set $J_{j}$ from (4): $J_{j}=\{i\in[n]\;:\;\mathbf{A}_{ji}\neq 0\}$.

Let $\hat{S}$ be an arbitrary sampling and $v=(v_{1},\dots,v_{n})$ be defined by:

Let $j\in[m]$ and $\mathbf{A}_{j:}$ denote the $j$th row vector of matrix $\mathbf{A}$. By the definition of $J_{j}$,

Thus, $\mathbf{P}(\hat{S})\circ(\mathbf{A}_{j:}^{\top}\mathbf{A}_{j:})=\mathbf{P}(\hat{S})\circ\mathbf{P}(\hat{E}_{J_{j}})\circ(\mathbf{A}_{j:}^{\top}\mathbf{A}_{j:})=\mathbf{P}(J_{j}\cap\hat{S})\circ(\mathbf{A}_{j:}^{\top}\mathbf{A}_{j:})$. We now apply Lemma 5.2 to the sampling $J_{j}\cap\hat{S}$ and the matrix $\mathbf{A}_{j:}$ and obtain:

where $p=(p_{1},\dots,p_{n})$ is the vector of probability defined in (1) and $v=(v_{1},\dots,v_{n})$ is defined in (58). For completeness, let us show that the second inequality in (59) can be replaced by equality. Indeed, from (40) and the fact that $|J_{j}|=\|A_{j:}\|_{0}$, we obtain $\lambda^{\prime}(\mathbf{A}_{j:}^{\top}\mathbf{A}_{j:})=|J_{j}|$. Finally, using the upper bound in Theorem 4.1, we know that $\lambda^{\prime}(J_{j}\cap\hat{S})\leqslant|J_{j}|$. Hence, $\min\{\lambda^{\prime}(J_{j}\cap\hat{S}),\lambda^{\prime}(\mathbf{A}_{j:}^{\top}\mathbf{A}_{j:})\}=\lambda^{\prime}(J_{j}\cap\hat{S})$. ∎

The benefit of this approach is twofold: First, if the data matrix $\mathbf{A}$ is sparse, the sets $J_{j}$ have small cardinality, and from Proposition 4.1 (or other results in Section 4.3, depending on the sampling $\hat{S}$ used) we conclude that $\lambda^{\prime}(J_{j}\cap\hat{S})$ is small. Hence, the parameters $v_{i}$ obtained through (58) get better (i.e., smaller) with sparser data. Second, the formula for $v_{i}$ does not involve the need to compute an eigenvalue associated with the data matrix. On the other hand, instead of having to compute $\lambda^{\prime}(\hat{S})$ (which, as we have seen, is equal to $\tau$ if $|\hat{S}|=\tau$ with probability 1), we now need to compute the normalized largest eigenvalue of $m$ restrictions of $\hat{S}$, $\lambda^{\prime}(J_{j}\cap\hat{S})$ for all $j=1,2,\dots,m$. However, for this there is a good upper bound available through Proposition 4.1 for an arbitrary sampling, and refined bounds can be derived for specific samplings (for examples, see Section 4.3).

4 ESO without eigenvalues

In this section we illustrate the use of the techniques developed in the preceding sections to derive ESO inequalities, for selected samplings, which do not depend on any eigenvalues, and lead to easily computable ESO parameters $v=(v_{1},\dots,v_{n})$. The techniques can be used to derive similar ESO inequalities for other samplings as well.

Let $f$ satisfy Assumption 2.1 and let sets $J_{1},\dots,J_{m}$ be defined as in (4). Then $(f,\hat{S})\sim ESO(v)$ provided that the sampling $\hat{S}$ and vector $v$ are chosen in any of the following ways:

$\hat{S}$ is an arbitrary sampling such that $|\hat{S}|\leqslant\tau$ with probability 1, and

$\hat{S}$ is the $(c,\tau)$-distributed sampling and

where $\omega_{j}^{\prime}\stackrel{{\scriptstyle\text{def}}}{{=}}|\{l:\mathcal{P}_{l}\cap J_{j}\neq 0\}|$ for $j\in[m]$.

$\hat{S}$ is the $\tau$-nice sampling (for $\tau\geqslant 1$) and

$\hat{S}$ is a doubly uniform sampling (which is not nil) and

$\hat{S}$ is a serial sampling (i.e., a sampling for which $|\hat{S}|=1$ with probability 1) and $v=(v_{1},\dots,v_{n})$ is defined as in (64).

A direct consequence of Theorem 5.2 and Proposition 4.1.

A direct consequence of Theorem 5.2 and Proposition 4.2.

This is a special case of part (ii) for $c=1$.

A direct consequence of Theorem 5.2 and Proposition 4.4.

For a graph sampling it is clear that $|J_{j}\cap\hat{S}|\leqslant 1$ with probability 1 for all $j\in[m]$. The result then follows from Theorem 5.2.

A special case of (v). Indeed, a single vertex is an independent set of a graph.

Remarks: Note that part (i) of Proposition 5.1 is a strict improvement on (57). Also, this is strict improvement, both in the quality of the bound and in generality of the sampling, on the result in , which was proved for uniform samplings only and where the bound involved $\max_{j}|J_{j}|$ instead of $|J_{j}|$. Part (ii) should be compared with the results obtained in and part (iii) with those in .

Discussion

As stressed before, smaller parameter $v=(v_{1},\dots,v_{n})$ leads to better convergence result (see Table 1) but computing the smallest admissible $v$ would require too large computational effort. Nevertheless, using a cheaply computed parameter $v=(v_{1},\dots,v_{n})$ would lead to large iteration complexity and slow convergence. The trade-off between the preprocessing time for computing the parameter $v=(v_{1},\dots,v_{n})$ and the iteration complexity of the algorithm shall be discussed next.

For specific samplings such as $\tau$-nice sampling and $(c,\tau)$-distributed sampling, admissible $v$ can be computed using dedicated formulae 61 and 62, which appeared respectively in and . For arbitrary sampling $\hat{S}$, admissible parameter $v$ can be computed according to 57, 58 or 60, which are given for the first time. While 58 requires computing the largest eigenvalue for $m$ matrices of sizes $\{J_{1},\dots,J_{m}\}$, both 57 and 60 can be computed in at most two passes over the data. In return, 58 provides a smaller parameter $v$ which improves the iteration complexity.

For approximating $\lambda^{\prime}(J_{j}\cap\hat{S})$, one can apply power method on the positive semidefinite matrix $\mathbf{P}(J_{j}\cap\hat{S})$. The number of operations needed in one iteration of the power method is $|J_{j}|^{2}$ and if we apply $T$ iterations of power method Note that as the matrix $\mathbf{P}(J_{j}\cap\hat{S})$ is positive semidefinite, the power method always converges to the largest eigenvalue even if it is not a dominant eigenvalue. We defer the study on the convergence rate of power method for different matrices $\mathbf{P}(J_{j}\cap\hat{S})$ to a future work. , then the total number of operations needed for computing $v$ using 58 is

where the big $O$ notation hides constants independent of the data matrix $\mathbf{A}$.

Recall from Table 1 how the iteration complexity of different methods depends on the parameter $v=(v_{1},\dots,v_{n})$. Let us consider the strongly convex smooth objective function setup and assume that the random sampling $\hat{S}$ has cardinality $\tau$ with probability 1. Then the computational time of one epoch ($n$ iterations) is of the same order as $\tau$ passes over the data. Therefore, given a parameter $v=(v_{1},\dots,v_{n})$, the number of passes over the data is bounded by:

where $\epsilon$ is the target accuracy and $\lambda$ is the strong convexity parameter of the problem.

The comparison of the three formulae in terms of overall complexity is reported in Table 2, where the big $O$ notation hides constants independent of the data matrix $\mathbf{A}$. It is clear from the table that the trade-off between the preprocessing and the iteration complexity mainly depends on the proportion between $\tfrac{T\sum_{j=1}^{m}|J_{j}|^{2}}{\operatorname{nnz}(\mathbf{A})}$ and $\frac{1}{\lambda}\max_{i}\frac{v_{i}\tau}{p_{i}n}$. In Table 3 we report the actual computing time of $v$ using different formulae and the corresponding value of $\max_{i}\frac{v_{i}\tau}{p_{i}n}$, for two real data matrices w8a and dorothea. To facilitate the comparison we normalized the two data sets so that the diagonal elements of $\mathbf{A}^{\top}\mathbf{A}$ are all one. The samplings $\hat{S}$ that we used in the experiments are all product sampling (Definition 2.8) with respect to some random partition of the set $[n]$. The number of iterations $T$ for the power method is fixed to 10 and we multiply the obtained value by 1.01. Because of the comparable processing time, Formula 60 is clearly better than Formula 57. From Table 3 we also see that Formula 58 requires significant computational effort for computing $v$ comparing to the other two but also reduces the value of $\max_{i}\frac{v_{i}\tau}{p_{i}n}$ by order of magnitude in most of the regimes. Let us take the example of dorothea with $\tau=256$, then the overall number of passes over data is $O(1.52+\frac{256.91}{\lambda}\log(\frac{1}{\epsilon}))$ if $v$ is computed using Formula 60 and $O(8715.8+\frac{16.68}{\lambda}\log(\frac{1}{\epsilon}))$ if $v$ is computed using Formula 58. Hence for small enough strong convexity parameter $\lambda$, it is worth to spend more time in computing a good parameter $v$ using Formula 58, which will then be compensated by a smaller iteration complexity.

2 Optimal sampling

Proposition 5.1 should be understood in the context of complexity results for randomized coordinate descent, such as those in Table 1. For instance, in view of (60) for an arbitrary sampling $\hat{S}$ such that $|\hat{S}|\leqslant\tau$ with probability 1, the accelerated coordinate descent method developed in has complexity

Naturally, the bound improves if we use a specialized sampling, such as the $\tau$-nice sampling (since the constants $v_{i}$ become smaller).

Sometimes, one can find a sampling which minimizes the complexity bound. For instance, if we restrict our attention to serial samplings only (samplings picking a single coordinate at a time), then one can find probabilities $p_{1},\dots,p_{n}$, which uniquely define a sampling, minimizing the complexity bound:

where $w_{i}=\sum_{j}\mathbf{A}_{ji}^{2}$. Note that if the $i$th coordinate is optimal at the starting point (i.e., if $x_{i}^{0}=x_{i}^{*}$), then the prediction is to choose $p_{i}=0$ (i.e., to never update coordinate $i$) – this is what one would expect. Using the serial sampling defined by (66), the complexity (65) takes the form

While $\|d\|_{6}\leqslant\|d\|_{2}$ for all $d$, these quantities can be equal, in which case $C_{opt}$ is $n$ times better than $C_{unif}$.

Conclusion

We have conducted a systematic study of ESO inequalities for a large class of functions (those satisfying Assumption 2.1) and arbitrary samplings. These inequalities are crucial in the design and complexity analysis of randomized coordinate descent methods. This led us to study standard and normalized largest eigenvalue of the Hadamard product of the probability matrix associated with a sampling and a certain positive semidefinite matrix containing the data defining the function. Using our approach we have established new ESO results and also re-derived ESO results already established in the literature (in the case of uniform samplings) via different techniques. Our approach can be used to derive further bounds for specific samplings and can potentially be of interest outside the domain of randomized coordinate descent.

References

Appendix A Frequently used notation