On Optimal Probabilities in Stochastic Coordinate Descent Methods

We propose and analyze a new parallel coordinate descent method---`NSync---in which at each iteration a random subset of coordinates is updated, in parallel, allowing for the subsets to be chosen non-uniformly. We derive convergence rates under a strong convexity assumption, and comment on how to assign probabilities to the sets to optimize the bound. The complexity and practical performance of the method can outperform its uniform variant by an order of magnitude. Surprisingly, the strategy of updating a single randomly selected coordinate per iteration---with optimal probabilities---may require less iterations, both in theory and practice, than the strategy of updating all coordinates at every iteration.


Introduction
In this work we consider the optimization problem where φ is strongly convex and smooth. We propose a new algorithm, and call it 'NSync (Nonuniform SYNchronous Coordinate descent).

Algorithm 1 ('NSync)
Input: Initial point x 0 ∈ R n , subset probabilities {p S } and stepsize parameters w 1 , . . . , w n > 0 for k = 0, 1, 2, . . . do Select a random set of coordinatesŜ ⊆ {1, . . . , n} such that Prob(Ŝ = S) = p S Updated selected coordinates: In 'NSync, we first assign a probability p S ≥ 0 to every subset S of [n] := {1, . . . , n}, with S p S = 1, and pick stepsize parameters w i > 0, i = 1, 2, . . . , n. At every iteration, a random setŜ is generated, independently from previous iterations, following the law Prob(Ŝ = S) = p S , and then coordinates i ∈Ŝ are updated in parallel by moving in the direction of the negative partial derivative with stepsize 1/w i . The updates are synchronized: no processor/thread is allowed to proceed before all updates are applied, generating the new iterate x k+1 . We specifically study samplingsŜ which are non-uniform in the sense that p i := Prob(i ∈Ŝ) = S:i∈S p S is allowed to vary with i. By ∇ i φ(x) we mean ∇φ(x), e i , where e i ∈ R n is the i-th unit coordinate vector.
Literature. Serial stochastic coordinate descent methods were proposed and analyzed in [6,13,15,18], and more recently in various settings in [12,7,8,9,21,19,24,3]. Parallel methods were considered in [2,16,14], and more recently in [22,5,23,4,11,20,10,1]. A memory distributed method scaling to big data problems was recently developed in [17]. A nonuniform coordinate descent method updating a single coordinate at a time was proposed in [15], and one updating two coordinates at a time in [12]. To the best of our knowledge, 'NSync is the first nonuniform parallel coordinate descent method.

Analysis
Our analysis of 'NSync is based on two assumptions. The first assumption generalizes the ESO concept introduced in [16] and later used in [22,23,5,4,17] to nonuniform samplings. The second assumption requires that φ be strongly convex.
Notation: For x, y, u ∈ R n we write . , x n y n ) and u −1 := (1/u 1 , . . . , 1/u n ). For S ⊆ [n] and h ∈ R n , let h [S] := i∈S h i e i . Assumption 1 (Nonuniform ESO: Expected Separable Overapproximation). Assume p = (p 1 , . . . , p n ) T > 0 and that for some positive vector w ∈ R n and all x, h ∈ R n , Inequalities of type (2), in the uniform case (p i = p j for all i, j), were studied in [16,22,5,17]. Assumption 2 (Strong convexity). We assume that φ is γ-strongly convex with respect to the norm We can now establish a bound on the number of iterations sufficient for 'NSync to approximately solve (1) with high probability. Theorem 3. Let Assumptions 1 and 2 be satisfied. Choose If {x k } are the random iterates generated by 'NSync, then Moreover, we have the lower bound Λ ≥ ( i wi vi )/E[|Ŝ|]. Proof. We first claim that φ is µ-strongly convex with respect to the norm · w•p −1 , i.e., where µ := γ/Λ. Indeed, this follows by comparing (3) and (6) in the light of (4). Let x * be such Let , and utilizing Assumption 1, we get Taking expectations in the last inequality and rearranging the terms, we obtain Using this, Markov inequality, and the definition of K, we finally get Let us now establish the last claim. First, note that (see [16,Sec 3.2] for more results of this type), Letting where the last equality follows since optimal p i is proportional to Theorem 3 is generic in the sense that we do not say when Assumptions 1 and 2 are satisfied, how should one go about to choose the stepsizes w and probabilities {p S }. In the next section we address these issues. On the other hand, this abstract setting allowed us to write a brief complexity proof.

Nonuniform samplings and ESO
Consider now problem (1) with φ of the form where v > 0. Note that Assumption 2 is satisfied. We further make the following two assumptions. Assumption 4 (Smoothness). f has Lipschitz gradient with respect to the coordinates, with positive constants L 1 , . . . , L n . That is, where J is a finite collection of nonempty subsets of [n] and f J are differentiable convex functions such that f J depends on coordinates i ∈ J only. Let ω := max J |J|. We say that f is separable of degree ω.
Uniform parallel coordinate descent methods for regularized problems with f of the above structure were analyzed in [16].
whence ω is the maximum # of nonzeros in a row of A.
Nonuniform sampling. Instead of considering the general case of arbitrary p S assigned to all subsets of [n], here we consider a special kind of sampling having two advantages: i) sets can be generated easily, ii) it leads to larger stepsizes 1/w i and hence improved convergence rate. Fix τ ∈ [n] and c ≥ 1 and let S 1 , . . . , S c be a collection of (possibly overlapping) subsets of [n] such that |S j | ≥ τ for all i and ∪ c j=1 S j = [n]. Moreover, let q = (q 1 , . . . , q c ) > 0 be a probability vector. LetŜ j be τ -nice sampling from S j ; that is,Ŝ j picks subsets of S j having cardinality τ , uniformly at random. We assume these samplings are independent. Now,Ŝ is generated as follows. We first pick j ∈ {1, . . . , c} with probability q j , and then drawŜ j . Note that we do not need to compute the quantities p S , S ⊆ [n], to execute 'NSync. In fact, it is much easier to implement the sampling via the two-tier procedure explained above. SamplingŜ is a nonuniform variant of the τ -nice sampling studied in [16], which here arises as a special case for c = 1. Note that where δ ij = 1 if i ∈ S j , and 0 otherwise.
Theorem 7. Let Assumptions 4 and 5 be satisfied, and letŜ be the sampling described above. Then Assumption 1 is satisfied with p given by (12) and any w = (w 1 , . . . , w n ) T for which where ω j := max J∈J |J ∩ S j | ≤ ω.
Proof. Since f is separable of degree ω, so is φ (because 1 2 x 2 v is separable). Now, where the last inequality follows from the ESO for τ -nice samplings established in [16,Theorem 15]. The claim now follows by comparing the above expression and (2).

Optimal probabilities
Observe that formula (13) can be used to design a sampling (characterized by the sets S j and probabilities q j ) that minimizes Λ, which in view of Theorem 3 optimizes the convergence rate of the method.
Serial setting. Consider the serial version of 'NSync (Prob(|Ŝ| = 1) = 1). We can model this via c = n, with S i = {i} and p i = q i for all i ∈ [n]. In this case, using (12) and (13), we get w i = w * i = L i + v i . Minimizing Λ in (4) over the probability vector p gives the optimal probabilities (we refer to this as the optimal serial method) and optimal complexity respectively. Note that the uniform sampling, p i = 1/n for all i, leads to Λ U S := n + n max j L j /v j (we call this the uniform serial method), which can be much larger than Λ Optimal serial method can be faster than the fully parallel method. To model the fully parallel setting (i.e., the variant of 'NSync updating all coordinates at every iteration), we can set c = 1 and τ = n, which yields Λ F P = ω + ω max j L j /v j . Since ω ≤ n, it is clear that Λ U S ≥ Λ F P . However, for large enough ω it will be the case that Λ F P ≥ Λ OS , implying, surprisingly, that the optimal serial method can be faster than the fully parallel method.
Parallel setting. Fix τ and sets S j , j = 1, 2, . . . , c, and define θ := max j 1 + i , so we are fine). From (4), (12) and (13) we see that the complexity of 'NSync is determined by The probability vector q minimizing this quantity can be computed by solving a linear program with c+1 variables (q 1 , . . . , q c , α), 2n linear inequality constraints and a single linear equality constraint:

Experiments
We now conduct 2 preliminary small scale experiments to illustrate the theory; the results are depicted below. All experiments are with problems of the form (11) with f chosen as in Example 6. In the left plot we chose A ∈ R 2×30 , γ = 1, v 1 = 0.05, v i = 1 for i = 1 and L i = 1 for all i. We compare the US method (p i = 1/n, blue) with the OS method (p i given by (16), red). The dashed lines show 95% confidence intervals (we run the methods 100 times, the line in the middle is the average behavior). While OS can be faster, it is sensitive to over/under-estimation of the constants L i , v i . In the right plot we show that a nonuniform serial (NS) method can be faster than the fully parallel (FP) variant (we have chosen m = 8, n = 10 and 3 values of ω). On the horizontal axis we display the number of epochs, where 1 epoch corresponds to updating n coordinates (for FP this is a single iteration, whereas for NS it corresponds to n iterations).