Stochastic Dual Coordinate Ascent with Adaptive Probabilities
Dominik Csiba, Zheng Qu, Peter Richtárik
Introduction
Empirical Loss Minimization. In this paper we consider the regularized empirical risk minimization problem:
We assume throughout that the loss functions are -smooth for some . That is, we assume they are differentiable and have Lipschitz derivative with Lipschitz constant :
for all , and .
The ERM problem (1) has received considerable attention in recent years due to its widespread usage in supervised statistical learning (Shalev-Shwartz & Zhang, 2013b). Often, the number of samples is very large and it is important to design algorithms that would be efficient in this regime.
Modern stochastic algorithms for ERM. Several highly efficient methods for solving the ERM problem were proposed and analyzed recently. These include primal methods such as SAG (Schmidt et al., 2013), SVRG (Johnson & Zhang, 2013), S2GD (Konečný & Richtárik, 2014), SAGA (Defazio et al., 2014), mS2GD (Konečný et al., 2014a) and MISO (Mairal, 2014). Importance sampling was considered in ProxSVRG (Xiao & Zhang, 2014) and S2CD (Konečný et al., 2014b).
Stochastic Dual Coordinate Ascent. One of the most successful methods in this category is stochastic dual coordinate ascent (SDCA), which operates on the dual of the ERM problem (1):
where functions and are defined by
SDCA in each iteration randomly selects a dual variable , and performs its update, usually via closed-form formula – this strategy is know as randomized coordinate descent. Methods based on updating randomly selected dual variables enjoy, in our setting, a linear convergence rate (Shalev-Shwartz & Zhang, 2013b, 2012; Takáč et al., 2013; Shalev-Shwartz & Zhang, 2013a; Zhao & Zhang, 2014; Qu et al., 2014). These methods have attracted considerable attention in the past few years, and include SCD (Shalev-Shwartz & Tewari, 2011), RCDM (Nesterov, 2012), UCDC (Richtárik & Takáč, 2014), ICD (Tappenden et al., 2013), PCDM (Richtárik & Takáč, 2012), SPCDM (Fercoq & Richtárik, 2013), SPDC (Zhang & Xiao, 2014), APCG (Lin et al., 2014), RCD (Necoara & Patrascu, 2014), APPROX (Fercoq & Richtárik, 2013), QUARTZ (Qu et al., 2014) and ALPHA (Qu & Richtárik, 2014). Recent advances on mini-batch and distributed variants can be found in (Liu & Wright, 2014), (Zhao et al., 2014b), (Richtárik & Takáč, 2013a), (Fercoq et al., 2014), (Trofimov & Genkin, 2014), (Jaggi et al., 2014), (Mareček et al., 2014) and (Mahajan et al., 2014). Other related work includes (Nemirovski et al., 2009; Duchi et al., 2011; Agarwal & Bottou, 2014; Zhao et al., 2014a; Fountoulakis & Tappenden, 2014; Tappenden et al., 2014). We also point to (Wright, 2014) for a review on coordinate descent algorithms.
Selection Probabilities. Naturally, both the theoretical convergence rate and practical performance of randomized coordinate descent methods depends on the probability distribution governing the choice of individual coordinates. While most existing work assumes uniform distribution, it was shown by Richtárik & Takáč (2014); Necoara et al. (2012); Zhao & Zhang (2014) that coordinate descent works for an arbitrary fixed probability distribution over individual coordinates and even subsets of coordinates (Richtárik & Takáč, 2013b; Qu et al., 2014; Qu & Richtárik, 2014, 2014). In all of these works the theory allows the computation of a fixed probability distribution, known as importance sampling, which optimizes the complexity bounds. However, such a distribution often depends on unknown quantities, such as the distances of the individual variables from their optimal values (Richtárik & Takáč, 2014; Qu & Richtárik, 2014). In some cases, such as for smooth strongly convex functions or in the primal-dual setup we consider here, the probabilities forming an importance sampling can be explicitly computed (Richtárik & Takáč, 2013b; Zhao & Zhang, 2014; Qu et al., 2014; Qu & Richtárik, 2014, 2014). Typically, the theoretical influence of using the importance sampling is in the replacement of the maximum of certain data-dependent quantities in the complexity bound by the average.
Adaptivity. Despite the striking developments in the field, there is virtually no literature on methods using an adaptive choice of the probabilities. We are aware of a few pieces of work; but all resort to heuristics unsupported by theory (Glasmachers & Dogan, 2013; Lukasewitz, 2013; Schaul et al., 2013; Banks-Watson, 2012; Loshchilov et al., 2011), which unfortunately also means that the methods are sometimes effective, and sometimes not. We observe that in the primal-dual framework we consider, each dual variable can be equipped with a natural measure of progress which we call “dual residue”. We propose that the selection probabilities be constructed based on these quantities.
Outline: In Section 2 we summarize the contributions of our work. In Section 3 we describe our first, theoretical methods (Algorithm 1) and describe the intuition behind it. In Section 4 we provide convergence analysis. In Section 5 we introduce Algorithm 2: an variant of Algorithm 1 containing heuristic elements which make it efficiently implementable. We conclude with numerical experiments in Section 6. Technical proofs and additional numerical experiments can be found in the appendix.
Contributions
We now briefly highlight the main contributions of this work.
Two algorithms with adaptive probabilities. We propose two new stochastic dual ascent algorithms: AdaSDCA (Algorithm 1) and AdaSDCA+ (Algorithm 2) for solving (1) and its dual problem (2). The novelty of our algorithms is in adaptive choice of the probability distribution over the dual coordinates.
Complexity analysis. We provide a convergence rate analysis for the first method, showing that AdaSDCA enjoys better rate than the best known rate for SDCA with a fixed sampling (Zhao & Zhang, 2014; Qu et al., 2014). The probabilities are proportional to a certain measure of dual suboptimality associated with each variable.
Practical method. AdaSDCA requires the same computational effort per iteration as the batch gradient algorithm. To solve this issue, we propose AdaSDCA+ (Algorithm 2): an efficient heuristic variant of the AdaSDCA. The computational effort of the heuristic method in a single iteration is low, which makes it very competitive with methods based on importance sampling, such as IProx-SDCA (Zhao & Zhang, 2014). We support this with computational experiments in Section 6.
Outline: In Section 2 we summarize the contributions of our work. In Section 3 we describe our first, theoretical methods (AdaSDCA) and describe the intuition behind it. In Section 4 we provide convergence analysis. In Section 5 we introduce AdaSDCA+: a variant of AdaSDCA containing heuristic elements which make it efficiently implementable. We conclude with numerical experiments in Section 6. Technical proofs and additional numerical experiments can be found in the appendix.
The Algorithm: AdaSDCA
where is the -by- matrix with columns .
Note, that if and only if satisfies (5). This motivates the design of AdaSDCA (Algorithm 1) as follows: whenever is large, the th dual coordinate is suboptimal and hence should be updated more often.
Alternatively, is coherent with if for
AdaSDCA is a stochastic dual coordinate ascent method, with an adaptive probability vector , which could potentially change at every iteration . The primal and dual update rules are exactly the same as in standard SDCA (Shalev-Shwartz & Zhang, 2013b), which instead uses uniform sampling probability at every iteration and does not require the computation of the dual residue .
Consider the AdaSDCA algorithm during iteration and assume that is coherent with . Then
Lemma 3 is proved similarly to Lemma 2 in (Zhao & Zhang, 2014), but in a slightly more general setting. For completeness, we provide the proof in the appendix. ∎
We also need to make sure that in order to apply Lemma 3. A “good” adaptive probability should then be the solution of the following optimization problem:
A feasible solution to (11) is the importance sampling (also known as optimal serial sampling) defined by:
which was proposed in (Zhao & Zhang, 2014) to obtain proximal stochastic dual coordinate ascent method with importance sampling (IProx-SDCA). The same optimal probability vector was also deduced, via different means and in a more general setting in (Qu et al., 2014). Note that in this special case, since is independent of the residue , the computation of is unnecessary and hence the complexity of each iteration does not scale up with .
It seems difficult to identify other feasible solutions to program (11) apart from , not to mention solve it exactly. However, by relaxing the constraint , we obtain an explicit optimal solution.
The optimal solution of
The suggestion made by (14) is clear: we should update more often those dual coordinates which have large absolute dual residue and/or large Lipschitz constant .
If we let and , the constraint (9) may not be sastified, in which case (8) does not necessarily hold. However, as shown by the next lemma, the constraint (9) is not required for obtaining (8) when all the functions are quadratic.
Suppose that all are quadratic. Let . If , then (8) holds for any .
Convergence results
In this section we present our theoretical complexity results for AdaSDCA. The main results are formulated in Theorem 7, covering the general case, and in Theorem 11 in the special case when are all quadratic.
We derive the convergence result from Lemma 3.
Let . If and , then
This follows directly from Lemma 3 and the fact that the right-hand side of (8) equals 0 when . ∎
Consider AdaSDCA. If at each iteration , and , then
By plugging the last bound into (17) we get the bound on the primal dual error:
As mentioned in Section 3, by letting every sampling probability be the importance sampling (optimal serial sampling) defined in (12), AdaSDCA reduces to IProx-SDCA proposed in (Zhao & Zhang, 2014). The convergence theory established for IProx-SDCA in (Zhao & Zhang, 2014), which can also be derived as a direct corollary of our Theorem 7, is stated as follows.
Consider AdaSDCA with defined in (12) for all . Then
The next corollary suggests that a better convergence rate than IProx-SDCA can be achieved by using properly chosen adaptive sampling probability.
However, solving (11) requires large computational effort, because of the dimension and the non-convex structure of the program. We show in the next section that when all the loss functions are quadratic, then we can get better convergence rate in theory than IProx-SDCA by using the optimal solution of (13).
2 Quadratic loss functions
The main difficulty of solving (11) comes from the inequality constraint, which originates from (9). In this section we mainly show that the constraint (9) can be released if all are quadratic.
Suppose that all are quadratic. Let . If , then
This is a direct consequence of Lemma 5 and the fact that the right-hand side of (8) equals 0 when . ∎
Suppose that all are quadratic. Consider AdaSDCA. If at each iteration , , then (15) holds for all .
We only need to apply Proposition 10. The rest of the proof is the same as in Theorem 7. ∎
Efficient heuristic variant
Corollary 9 and 12 suggest how to choose adaptive sampling probability in AdaSDCA which yields a theoretical convergence rate at least as good as IProx-SDCA (Zhao & Zhang, 2014). However, there are two main implementation issues of AdaSDCA:
The update of the dual residue at each iteration costs where is the number of nonzero elements of the matrix ;
We do not know how to compute the optimal solution of (11).
In this section, we propose a heuristic variant of AdaSDCA, which avoids the above two issues while staying close to the ’good’ adaptive sampling distribution.
AdaSDCA+ has the same structure as AdaSDCA with a few important differences.
Epochs AdaSDCA+ is divided into epochs of length . At the beginning of every epoch, sampling probabilities are computed according to one of two options. During each epoch the probabilities are cheaply updated at the end of every iteration to approximate the adaptive model. The intuition behind is as follows. After is sampled and the dual coordinate is updated, the residue naturally decreases. We then decrease also the probability that is chosen in the next iteration, by setting to be proportional to . By doing this we avoid the computation of at each iteration (issue 1) which costs as much as the full gradient algorithm, while following closely the changes of the dual residue . We reset the adaptive sampling probability after every epoch of length .
Parameter The setting of parameter in AdaSDCA+ directly affects the performance of the algorithm. If is too large, the probability of sampling the same coordinate twice during an epoch will be very small. This will result in a random permutation through all coordinates every epoch. On the other hand, for too small the coordinates having larger probabilities at the beginning of an epoch could be sampled more often than it should, even after their corresponding dual residues become sufficiently small. We don’t have a definitive rule on the choice of and we leave this to future work. Experiments with different choices of can be found in Section 6.
Option I & Option II At the beginning of each epoch, one can choose between two options for resetting the sampling probability. Option I corresponds to the optimal solution of (13), given by the closed form (14). Option II is the optimal serial sampling probability (12), the same as the one used in IProx-SDCA (Zhao & Zhang, 2014). However, AdaSDCA+ differs significantly with IProx-SDCA since we also update iteratively the sampling probability, which as we show through numerical experiments yields a faster convergence than IProx-SDCA.
2 Computational cost
Sampling and probability update During the algorithm we sample from non-uniform probability distribution , which changes at each iteration. This process can be done efficiently using the Random Counters algorithm introduced in Section 6.2 of (Nesterov, 2012), which takes operations to create the probability tree and operations to sample from the distribution or change one of the probabilities.
Total computational cost We can compute the computational cost of one epoch. At the beginning of an epoch, we need operations to calculate the dual residue . Then we create a probability tree using operations. At each iteration we need operations to sample a coordinate, operations to calculate the update to and a further operations to update the probability tree. As a result an epoch needs operations. For comparison purpose we list in Table 1 the one epoch computational cost of comparable algorithms.
Numerical Experiments
In this section we present results of numerical experiments.
In both cases we use -regularizer, i.e.,
Quadratic loss functions appear usually in regression problems, and smoothed Hinge loss can be found in linear support vector machine (SVM) problems (Shalev-Shwartz & Zhang, 2013a).
2 Numerical results
We used 5 different datasets: w8a, dorothea, mushrooms, cov1 and ijcnn1 (see Table 2).
In all our experiments we used and .
AdaSDCA The results of the theory developed in Section 4 can be observed through Figure 1 to Figure 4. AdaSDCA needs the least amount of iterations to converge, confirming the theoretical result.
AdaSDCA+ V.S. others We can observe through Figure 15 to 24, that both options of AdaSDCA+ outperforms SDCA and IProx-SDCA, in terms of number of iterations, for quadratic loss functions and for smoothed Hinge loss functions. One can observe similar results in terms of time through Figure 5 to Figure 14.
Option I V.S. Option II Despite the fact that Option I is not theoretically supported for smoothed hinge loss, it still converges faster than Option II on every dataset and for every loss function. The biggest difference can be observed on Figure 13, where Option I converges to the machine precision in just 15 seconds.
Different choices of To show the impact of different choices of on the performance of AdaSDCA+, in Figures 25 to 33 we compare the results of the two options of AdaSDCA+ using different equal to , and . It is hard to draw a clear conclusion here because clearly the optimal shall depend on the dataset and the problem type.
References
Proofs
It can be easily checked that the following relations hold
where is the output sequence of Algorithm 1. Let and . For each , since is -smooth, is -strongly convex and thus for arbitrary ,
where the last equality follows from the definition of in Algorithm 1. Then by letting for some arbitrary we get:
By taking expectation with respect to we get:
Then for each and by plugging it into (23) we get:
Note that (13) is a standard constrained maximization problem, where everything independent of can be treated as a constant. We define the Lagrangian
and get the following optimality conditions: