We analyze fast stochastic algorithms for the following optimization problem:
where is a Riemannian manifold equipped with the metric , and
is a random variable. We assume that for any, the function is geodesically -smooth (see Section 2
). This class of functions includes as special cases important problems such as principal component analysis (PCA), independent component analysis (ICA), dictionary learning, mixture modeling, among others. Moreover, the finite-sum problem () is a special case in which finite number of component functions are chosen uniformly at random (e.g., in Empirical Risk Minimization).
When solving problems with parameters constrained to lie on a manifold, a naive approach is to alternate between optimizing the cost in a suitable ambient Euclidean space and “projecting” onto the manifold. For example, two well-known methods to compute the leading eigenvector of symmetric matrices, power iteration and Oja’s algorithm(Oja, 1992), are in essence projected gradient and projected stochastic gradient algorithms. For certain manifolds (e.g., positive definite matrices), projections can be quite expensive to compute and possibly the Euclidean approach may have poor numerical conditioning (Yuan et al., 2017).
An effective alternative is to use Riemannian optimization, which directly operates on the manifold in question. This mode of operation allows Riemannian optimization to view the constrained optimization problem (1) as an unconstrained problem on a manifold, and thus, to be “projection-free.” More important is the conceptual viewpoint: by casting the problem in a Riemannian framework, one can discover insights into problem geometry that can translate into not only more precise mathematical analysis but also more efficient optimization algorithms.
The Euclidean version of (1) where and
is the Euclidean inner-product has been the subject of intense algorithmic development in machine learning and optimization, starting with the classical work ofRobbins and Monro (1951). However, both batch and stochastic gradient methods suffer from high computation load. For solving finite sum problems with components, the full-gradient method requires derivatives at each step; the stochastic method requires only one derivative, but at the expense of slower convergence to an
-accurate solution. These issues have motivated much of the progress on faster stochastic optimization in vector spaces by using variance reduction(Schmidt et al., 2013; Johnson and Zhang, 2013; Defazio et al., 2014; Konečnỳ and Richtárik, 2013). Along with many recent works (see related work), these algorithms achieve faster convergence than the original gradient descent algorithms in multiple settings.
Riemannian counterparts of batch and stochastic optimization algorithms have witnessed growing interest recently. Zhang and Sra (2016) present the first global complexity analysis of batch and stochastic gradient methods for geodesically convex functions. Later work (Zhang et al., 2016; Kasai et al., 2016; Sato et al., 2017) improves the convergence rate for finite-sum problems by using variance reduction techniques. In this paper, we develop this line of work further, and improve rates based on a more careful control of variance analyzed in (Fang et al., 2018; Nguyen et al., 2017). Another important aspect of our work is that by pursuing a slightly different analysis, we are able to remove the assumption that all iterates remain in a compact subset of the Riemannian manifold. Such an assumption was crucial to most previous Riemannian methods, but was not always fully justified.
We summarize the key contributions of this paper below.
We introduce R-Spider, a Riemannian variance reduced stochastic gradient method based on the recent SPIDER algorithm (Fang et al., 2018). We analyze R-Spider
for optimizing geodesically smooth stochastic nonconvex functions. To our knowledge, we obtain the first rate faster than Riemannian stochastic gradient descent for general nonconvex stochastic Riemannian optimization.
We specialize R-Spider to (Riemannian) nonconvex finite-sum problems. Our rate improves the best known rates and match the lower bound as in the Euclidean case.
We propose two variations of R-Spider for geodesically strongly convex problems and for Riemannian gradient dominated costs. For these settings, we achieve the best known rates in terms of number of samples and the condition number .
Importantly, we provide convergence guarantees that are independent of the Riemannian manifold’s diameter and its sectional curvature. This contribution is important in two main aspects. First, the best known theoretical upper bounds are improved. Second, the algorithm no longer assumes bounded diameter of the Riemannian manifold, which helps lift the assumption crucial for previous work that required all the iterates generated by the algorithm to remain in a compact set.
We briefly summarize the rates obtained in Table 1.
1.1 Related Work
Variance reduction in stochastic optimization.
Variance reduction techniques, such as control variates, are widely used in Monte Carlo simulations (Rubinstein and Kroese, 2011). In linear spaces, variance reduced methods for solving finite-sum problems have recently witnessed a huge surge of interest (e.g. Schmidt et al., 2013; Johnson and Zhang, 2013; Defazio et al., 2014; Bach and Moulines, 2013; Konečnỳ and Richtárik, 2013; Xiao and Zhang, 2014; Gong and Ye, 2014). They have been shown to accelerate finite sum optimization for strongly convex objectives and convex objectives. Later work by Lin et al. (2015); Allen-Zhu (2017a) further accelerates the rates in convex problems using techniques similar to Nesterov’s acceleration method (Nesterov, 2013). For nonconvex problems, Reddi et al. (2016); Allen-Zhu (2017b); Lei et al. (2017); Fang et al. (2018); Nguyen et al. (2017) also achieved faster rate than the vanilla (stochastic) gradient descent method in both finite sum settings and stochastic settings. Our analysis is inspired mainly by Fang et al. (2018); Zhang et al. (2016). The analysis can also be applied to (Wang et al., 2018) and achieve matching rate assuming access to proximal oracle on Riemannian manifold.
Earlier references can be found in Udriste (1994); Absil et al. (2009), where analysis is limited to asymptotic convergence (except (Udriste, 1994, Theorem 4.2)). Stochastic Riemannian optimization has been previously considered in Bonnabel (2013); Liu et al. (2004), though with only asymptotic convergence analysis, and without any rates. Many applications of Riemannian optimization are known, including matrix factorization on fixed-rank manifold (Vandereycken, 2013; Tan et al., 2014), dictionary learning (Cherian and Sra, 2015; Sun et al., 2015), optimization under orthogonality constraints (Edelman et al., 1998; Moakher, 2002)
, covariance estimation(Wiesel, 2012), learning elliptical distributions (Zhang et al., 2013; Sra and Hosseini, 2013), Poincaré embeddings (Nickel and Kiela, 2017)2015). Zhang and Sra (2016) provide the first global complexity analysis for first-order Riemannian algorithms, but their analysis is restricted to geodesically convex problems with full or stochastic gradients. Boumal et al. (2016) analyzed iteration complexity of Riemannian trust-region methods, whereas Bento et al. (2017) studied non-asymptotic convergence of Riemannian gradient, subgradient and proximal point methods. Tripuraneni et al. (2018); Zhang and Sra (2018) analyzed aspects other than variance reduction to accelerate the convergence of first order optimization methods on Riemannian manifolds. Zhang et al. (2016); Sato et al. (2017) analyzed variance reduction techniques on Riemannian manifolds, and their rate has remain best-known up to our knowledge. In this work, we improve upon their results. Zhou et al. (2018) worked on the same problem in parallel and achieved the same rate. The difference between this work and Zhou et al. (2018) is that our algorithm uses a constant step size and adaptive sample size. This enables us to bound instead of . Hence our result is slightly stronger and further simplifies later proof for the convergence of gradient dominated functions.
Before formally discussing Riemannian optimization, let us recall some foundational concepts of Riemannian geometry. For a thorough review one can refer to any classic text, e.g.,(Petersen, 2006).
A Riemannian manifold is a real smooth manifold equipped with a Riemannain metric . The metric induces an inner product structure in each tangent space associated with every . We denote the inner product of as ; and the norm of is defined as . The angle between is defined as . A geodesic is a constant speed curve that is locally distance minimizing. An exponential map maps in to on , such that there is a geodesic with and . If between any two points in there is a unique geodesic, the exponential map has an inverse and the geodesic is the unique shortest path with the geodesic distance between .
Parallel transport maps a vector to , while preserving norm, and roughly speaking, “direction,” analogous to translation in . A tangent vector of a geodesic remains tangent if parallel transported along . Parallel transport preserves inner products.
We now define some key terms. A set is called geodesically convex if for any , there is a geodesic with and for . Throughout the paper, we assume that the function in (1) is defined on a Riemannian manifold .
In the following we do not explicitly write Riemannian metric or the index of tangent space to simplify notation, as they should be obvious from the context: inner product of is defined as ; norm of is defined as .
Based on the above notations, we define the following properties of the function in (1).
Definition 1 (Strong convexity).
A function is said to be geodesically -strongly convex if for any ,
Definition 2 (Smoothness).
A differentiable function is said to be geodesically -smooth if its gradient is -Lipschitz, i.e. for any ,
where is the parallel transport from to .
Observe that compared to the Euclidean setup, the above definition requires a parallel transport operation to “transport” to . It can be proved that if is -smooth, then for any ,
Definition 3 (PL inequality).
is -gradient dominated if is a global minimizer of and for every
As in the Euclidean case, gradient dominated is implied by strongly convex.
An Incremental First-order Oracle (IFO) (Agarwal and Bottou, 2015) in (1) takes in a point , and generates a random sample . The oracle then returns a pair . In finite-sum setting, takes values in and each random sample corresponds to one of component functions. We measure non-asymptotic complexity in terms of IFO calls.
3 Riemannian SPIDER
In this section we introduce the R-SPIDER algorithm. In particular, we propose one variant for nonconvex problems, and two for gradient-dominated problems. Each variation aims to optimize a particular dependency on function parameters. Our proposed algorithm differs from the Euclidean SPIDER (Fang et al., 2018) in two key aspects: the variance reduction step uses parallel transport to combine gradients from different tangent spaces; and the exponential map is used (instead of the update ).
We would like to point out that if retractions instead of exponential maps are used in the proposed algorithms, our analysis will still hold if such that , , where .
3.1 General smooth nonconvex functions
Our proposed algorithm for solving nonconvex Riemannian optimization problems is shown in Algorithm 1. denotes the unbiased gradient estimator obtained by averaging samples. We first analyze the global complexity for solving nonconvex stochastic Riemannian optimization problems. In particular, we make the following assumptions
Assumption 1 (Smoothness).
For any fixed is geodesically smooth in .
Assumption 2 (Bounded objective).
Function is bounded below. Define where
Assumption 3 (Bounded variance).
Under the above assumptions, we make the following choice of parameters for running Algorithm 1.
We now state the following theorem for optimizing stochastic nonconvex functions.
Theorem 1 (Stochastic objective).
See Appendix A. The gist of our proof is to show that with sufficiently small variance of the gradient estimate, the algorithm either substantially decreases the objective function every iterations, or terminates because the gradient norm is small. ∎
Then we study the nonconvex problem under the finite-sum setting. In this setting, we assume . Hence we can write
We further make the following choice of parameters:
Then we have the following gurrantee.
Theorem 2 (Finite-sum objective).
The proof of the two theorems in this section follows by carefully applying the variance reduction technique proposed in Fang et al. (2018) onto the Riemannian manifold using the tools in Zhang et al. (2016). Unlike the SVRG-like algorithms in Zhang et al. (2016); Kasai et al. (2016), we can avoid analyzing the term , where is the snapshot point. Consequently, we do not need to resort to the trigonometric distance bound (see Zhang and Sra (2016)) and the convergence rate doesn’t depend on the sectional curvature bound.
Further, the convergence rates in both cases match their Euclidean counterparts. Remarkably, the rate under the finite-sum setting meets the lower bound as proved by Fang et al. (2018).
3.2 Gradient-dominated functions
In this section, we study the finite-sum problems with the following assumption.
is gradient dominated.
We denote as the condition number. To solve such problems, we propose two algorithms. The first algorithm is shown in Algorithm 2. It follows the same idea as in Zhang et al. (2016); Reddi et al. (2016). We have the following theorem on its convergence rate.
The second algorithm, shown in Algorithm 3, aims to achieve better complexity dependency on . With the following choice of parameters
we can make the following statement.
See Appendix D. The algorithm adaptively choose sample sizes based on the distance of the last update. The expected number of samples queried can be bounded by total sum of squared distances, which can further be bounded by the change in the objective value. ∎
In summary, Algorithm 2 achieves IFO complexity , while Algorithm 3 achieves sample complexity . It is unclear to us whether there exists an algorithm that performs uniformly better than both of the proposed algorithms. Further, we wish to point out that, as strong convexity implies gradient dominance, the convergence rates for the above algorithms also apply to strongly g-convex functions.
We introduce Riemannian SPIDER algorithms, a fast variance reduced stochastic gradient algorithm for Riemannian optimization. We analyzed the convergence rates of these algorithms for general smooth geodesically nonconvex functions under both finite-sum and stochastic settings, as well as for gradient dominated functions under the finite-sum setting. We showed that these algorithms improved the best known IFO complexity. We also removed the iteration complexity dependency on the curvature of the manifold.
There are a few open problems. First, the original SPIDER algorithm in Fang et al. (2018) and Algorithm 1 require very small stepsize. In practice, this usually results in very slow convergence rate. Even though the SPIDER-boost algorithm (Wang et al., 2018) and Algorithm 3 utilizes a constant large stepsize, the former one requires random termination of the algorithm, while the latter one requires very large sample size in each iteration. Therefore, none of these algorithms tend to perform well in practice if the implementation follows the theory exactly. Designing and testing practical algorithms with nice theoretical guarantees is left as future work.
Further, we approached the gradient-dominated functions with two different algorithms and got two different convergence rates. We suspect that it is not possible to achieve the best of both worlds at the same time. Proving such a lower bound or finding an algorithm that uniformly dominates both algorithms are both interesting research topics.
- Absil et al.  P.-A. Absil, R. Mahony, and R. Sepulchre. Optimization algorithms on matrix manifolds. Princeton University Press, 2009.
- Agarwal and Bottou  A. Agarwal and L. Bottou. A lower bound for the optimization of finite sums. In Proceedings of the 32nd International Conference on Machine Learning (ICML-15), pages 78–86, 2015.
Katyusha: The first direct acceleration of stochastic gradient
Proceedings of the 49th Annual ACM SIGACT Symposium on Theory of Computing, pages 1200–1205. ACM, 2017a.
- Allen-Zhu [2017b] Z. Allen-Zhu. Natasha 2: Faster non-convex optimization than sgd. arXiv preprint arXiv:1708.08694, 2017b.
- Bach and Moulines  F. Bach and E. Moulines. Non-strongly-convex smooth stochastic approximation with convergence rate o (1/n). In Advances in Neural Information Processing Systems, pages 773–781, 2013.
- Bento et al.  G. C. Bento, O. P. Ferreira, and J. G. Melo. Iteration-complexity of gradient, subgradient and proximal point methods on riemannian manifolds. Journal of Optimization Theory and Applications, 173(2):548–562, 2017.
- Bonnabel  S. Bonnabel. Stochastic gradient descent on Riemannian manifolds. Automatic Control, IEEE Transactions on, 58(9):2217–2229, 2013.
- Boumal et al.  N. Boumal, P.-A. Absil, and C. Cartis. Global rates of convergence for nonconvex optimization on manifolds. IMA Journal of Numerical Analysis, 2016.
- Cherian and Sra  A. Cherian and S. Sra. Riemannian dictionary learning and sparse coding for positive definite matrices. arXiv:1507.02772, 2015.
- Defazio et al.  A. Defazio, F. Bach, and S. Lacoste-Julien. Saga: A fast incremental gradient method with support for non-strongly convex composite objectives. In NIPS, pages 1646–1654, 2014.
- Edelman et al.  A. Edelman, T. A. Arias, and S. T. Smith. The geometry of algorithms with orthogonality constraints. SIAM journal on Matrix Analysis and Applications, 20(2):303–353, 1998.
- Fang et al.  C. Fang, C. J. Li, Z. Lin, and T. Zhang. Spider: Near-optimal non-convex optimization via stochastic path integrated differential estimator. arXiv preprint arXiv:1807.01695, 2018.
- Gong and Ye  P. Gong and J. Ye. Linear convergence of variance-reduced stochastic gradient without strong convexity. arXiv preprint arXiv:1406.1102, 2014.
- Hosseini and Sra  R. Hosseini and S. Sra. Matrix manifold optimization for Gaussian mixtures. In NIPS, 2015.
- Johnson and Zhang  R. Johnson and T. Zhang. Accelerating stochastic gradient descent using predictive variance reduction. In Advances in Neural Information Processing Systems, pages 315–323, 2013.
- Kasai et al.  H. Kasai, H. Sato, and B. Mishra. Riemannian stochastic variance reduced gradient on grassmann manifold. arXiv preprint arXiv:1605.07367, 2016.
- Konečnỳ and Richtárik  J. Konečnỳ and P. Richtárik. Semi-stochastic gradient descent methods. arXiv:1312.1666, 2013.
- Lei et al.  L. Lei, C. Ju, J. Chen, and M. I. Jordan. Non-convex finite-sum optimization via scsg methods. In Advances in Neural Information Processing Systems, pages 2348–2358, 2017.
- Lin et al.  H. Lin, J. Mairal, and Z. Harchaoui. A universal catalyst for first-order optimization. In Advances in Neural Information Processing Systems, pages 3384–3392, 2015.
- Liu et al.  X. Liu, A. Srivastava, and K. Gallivan. Optimal linear representations of images for object recognition. IEEE TPAMI, 26(5):662–666, 2004.
- Moakher  M. Moakher. Means and averaging in the group of rotations. SIAM journal on matrix analysis and applications, 24(1):1–16, 2002.
- Nesterov  Y. Nesterov. Introductory lectures on convex optimization: A basic course, volume 87. Springer Science & Business Media, 2013.
- Nguyen et al.  L. M. Nguyen, J. Liu, K. Scheinberg, and M. Takáč. Sarah: A novel method for machine learning problems using stochastic recursive gradient. arXiv preprint arXiv:1703.00102, 2017.
- Nickel and Kiela  M. Nickel and D. Kiela. Poincaré embeddings for learning hierarchical representations. In Advances in neural information processing systems, pages 6338–6347, 2017.
Principal components, minor components, and linear neural networks.Neural Networks, 5(6):927–935, 1992.
- Petersen  P. Petersen. Riemannian geometry, volume 171. Springer Science & Business Media, 2006.
- Reddi et al.  S. J. Reddi, A. Hefny, S. Sra, B. Póczós, and A. Smola. Stochastic variance reduction for nonconvex optimization. arXiv:1603.06160, 2016.
- Robbins and Monro  H. Robbins and S. Monro. A stochastic approximation method. Annals of Mathematical Statistics, 22:400–407, 1951.
- Rubinstein and Kroese  R. Y. Rubinstein and D. P. Kroese. Simulation and the Monte Carlo method, volume 707. John Wiley & Sons, 2011.
- Sato et al.  H. Sato, H. Kasai, and B. Mishra. Riemannian stochastic variance reduced gradient. arXiv preprint arXiv:1702.05594, 2017.
- Schmidt et al.  M. Schmidt, N. L. Roux, and F. Bach. Minimizing finite sums with the stochastic average gradient. arXiv:1309.2388, 2013.
- Sra and Hosseini  S. Sra and R. Hosseini. Geometric optimisation on positive definite matrices for elliptically contoured distributions. In Advances in Neural Information Processing Systems, pages 2562–2570, 2013.
- Sun et al.  J. Sun, Q. Qu, and J. Wright. Complete dictionary recovery over the sphere ii: Recovery by Riemannian trust-region method. arXiv:1511.04777, 2015.
- Tan et al.  M. Tan, I. W. Tsang, L. Wang, B. Vandereycken, and S. J. Pan. Riemannian pursuit for big matrix recovery. In International Conference on Machine Learning (ICML-14), pages 1539–1547, 2014.
- Tripuraneni et al.  N. Tripuraneni, N. Flammarion, F. Bach, and M. I. Jordan. Averaging stochastic gradient descent on Riemannian manifolds. arXiv preprint arXiv:1802.09128, 2018.
- Udriste  C. Udriste. Convex functions and optimization methods on Riemannian manifolds, volume 297. Springer Science & Business Media, 1994.
- Vandereycken  B. Vandereycken. Low-rank matrix completion by Riemannian optimization. SIAM Journal on Optimization, 23(2):1214–1236, 2013.
- Wang et al.  Z. Wang, K. Ji, Y. Zhou, Y. Liang, and V. Tarokh. Spiderboost: A class of faster variance-reduced algorithms for nonconvex optimization. arXiv preprint arXiv:1810.10690, 2018.
- Wiesel  A. Wiesel. Geodesic convexity and covariance estimation. IEEE Transactions on Signal Processing, 60(12):6182–6189, 2012.
- Xiao and Zhang  L. Xiao and T. Zhang. A proximal stochastic gradient method with progressive variance reduction. SIAM Journal on Optimization, 24(4):2057–2075, 2014.
- Yuan et al.  X. Yuan, W. Huang, P.-A. Absil, and K. A. Gallivan. A Riemannian quasi-Newton method for computing the Karcher mean of symmetric positive definite matrices. Florida State University, (FSU17-02), 2017.
- Zhang and Sra  H. Zhang and S. Sra. First-order methods for geodesically convex optimization. arXiv:1602.06053, 2016.
- Zhang and Sra  H. Zhang and S. Sra. Towards Riemannian accelerated gradient methods. arXiv preprint arXiv:1806.02812, 2018.
- Zhang et al.  H. Zhang, S. J. Reddi, and S. Sra. Riemannian SVRG: Fast stochastic optimization on Riemannian manifolds. In Advances in Neural Information Processing Systems, pages 4592–4600, 2016.
Zhang et al. 
T. Zhang, A. Wiesel, and M. S. Greco.
Multivariate generalized Gaussian distribution: Convexity and graphical models.Signal Processing, IEEE Transactions on, 61(16):4141–4148, 2013.
- Zhou et al.  P. Zhou, X.-T. Yuan, and J. Feng. Faster first-order methods for stochastic non-convex optimization on riemannian manifolds. arXiv preprint arXiv:1811.08109, 2018.
Appendix A Proof of Theorem 1
The proof of Theorem 1, Theorem 2 and Theorem 4 all follows from three steps: bound the variance of the gradient estimator; prove that function value decrease in expectation per iteration; bound the number of iterations with Assumption 2.
First we bound the variance of the estimator at each step.
Let be the sigma field generated by the random variable . Then forms a filtration. We can write the following equations
The first equality follows by the identities
The second equality follows by the fact that and
are unbiased estimators. Denotewhere is a sampled function from the distribution of as defined in (1). Note that . Hence we have
Substitue in (8) and we get that
where the first inequality follows by . The last inequality follows by the value of . Apply the bound recursively and denote , we get
Second, we show that the function value decreases.
By geodesically smoothness and (2), we have
The second inequality follows from the update rule of in Algorithm 1. The last inequality follows from . After taking expectation, we get
The inequality follows by .
Finally we are ready to prove the theorem.
Appendix B Proof of Theorem 2
The proof techniques are exactly the same as those in Section A, with the only changes being that and we use (6) rather than (4) as the parameters. We state two corresponding lemmas and leave out the details to avoid repetition.
Appendix C Proof of Theorem 3
By Theorem 2, we know that when the R-SPIDER-nonconvex algorithm terminates in iteration , it returns satisfying
By Assumption 4, we know that
By theorem 2, in iteration the R-SPIDER-nonconvex algorithm makes less than
IFO calls in expectation. Take expectation and substitute in the bound in (15), we get that the expected number of IFO calls is less than . ∎
Appendix D Proof of Theorem 4
As usual, we start with bounding the variance of the gradient estimator.
Then we prove that function value decreases in each epoch (iterations).